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ABSTRACT 


The indirect effects on the geoid computation due to the second 
method of Helmert's condensation were studied. When Helmert's anomalies 
are used in Stokes' equation, there are three types of corrections to the 
free-air geoid. The first correction, the indirect effect on geoid un- 
dulation due to the potential change in Helmert's reduction, had a max- 
imum value of 0.51 meters in the test area covering the United States. 

The second correction, the attraction change effect on geoid undulation, 
had a maximum value of 9.50 meters when the 10° cap was used in Stokes' 
equation. The last correction, the secondary indirect effect on geoid 
undulation, was found negligible in the test area. The corrections were 
applied to uncorrected free-air geoid undulations at 65 Doppler stations 
in the test area and cciiipared with the Doppler undulations. Based on the 
assumption that the Doppler coordinate system has a z-shift of 4 meters 
with respect to the geocenter, these comparisons showed that the cor- 
rections presented in this study yielded improved values of gravimetric 
undulations. 
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N geoid undulation 

Ag gravity anomaly 

h height above the geoid 

R mean radius of the earth 

y normal gravity at the point on the reference ellipsoid 

S (ip ) Stokes' function 

da element of solid angle 

W,V gravity potential and gravitational potential, respectively 

6W,6V gravity and gravitational potential change, respectively 

6N indirect effect due to the potential change 

AN attraction change effect on geoid undulation 

6 secondary indirect effect on gravity 

A,6A vertical attraction and attraction change, respectively 
C terrain correction 

k Newton’s gravitational constant 

p volume density 

k surface density 

\p angular radius from the evaluation point 

a linear radius from the evaluation point 

( )' quantity referred to the regular part of topography 

( )" quantity referred to the irregular part of topography 

( quantity referred to the topography before the condensation 

( )£or( ) s quantity referred to the condensed layer, i.e., after the condensation 

P point on the earth's surface (see figure 1) 

P Q point on the geoid (see figure 1) 

Q 0 point on the reference ellipsoid corresponding to point P 0 (see 

figure 1) 

C nm ,S nm fully normalized geopotential coefficients 
SBP spherical Bouguer plate (figure 4, p. 22) 

PCP plane circular plate (figure 5, p. 23) 

TDSBP Two Difference Spherical Bouguer Plate (p. 35) 

TDPCP Two Difference Plane Circular Plate (p. 35) 
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1. Introduction 


A fundamental equation used to compute geoid undulation, N is the 
Stokes' integral equation; 

N = ¥7 (/ Ag da (!) 

Equation (1) is valid for a reference ellipsoid which (i) has the same poten- 
tial as the geoid, (ii) encloses a mass that is numerically equal to the 
earth's mass, and (iii) has its center at the center of mass of the earth 
(Heiskanen and Moritz, 1967, p. 94). This implies that the gravity anomaly, 

Ag used in equation (1) refers to the geoid with no external 
masses. This important assumption necessitates the real earth to be regu- 
larized so that all masses outside the geoid are completely removed or shifted 
below the geoid. The process of removing or shifting the masses is called 
a gravity reduction. 

To properly obtain the gravity anomaly on the geoid, we must consider 
the effect of removing or shifting the masses outside the geoid on the value 
of gravity g obs at point P on the earth's surface. Then, after removing 
or shifting the masses, the gravity station is brought down along the vertical 
curve to point P 0 on the geoid (see Figure 1). In other words, the gravity 
force is transferred from point P to point P Q by a free air reduction. 

By subtracting the normal gravity y at the corresponding point, Q 0 , on 
the ellipsoid from the gravity at P 0 we then obtain the gravity anomaly 
on the geoid. 

By removing or shifting the masses outside the geoid in the gravity reduc- 
tion process the original potential of the earth is changed as well. Sup- 
pose that the geoid is defined to be' the equipotential surface having potential 
W 0 and the reduction process causes the change in potential at a point 
originally on the geoid by the amount of 5W . To be precise, SW is de- 
fined to be the gravity potential of the actual topographic masses minus the 
gravity potential of the masses after the reduction process. This means that 
the point originally on the geoid now has the potential W Q - 6W (the negative 
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sign is conventional, see below) and by definition it is not the point on 
the geoid anymore. We call the changed geoid as the cogeoid, Hence 
the surface computed by Stokes’ equation using the gravity anomaly reduced 
by the explained process so far is not the geoid but the cogeoid. The 
real geoid and the cogeoid are separated by a vertical distance 6N which 
is called an indirect effect of gravity reduction on geoid undulation. 
tSP! is simply computed by Brun's formula: 

«N = (2) 

where 6W is the gravity potential change at the geoid and y is the 
normal gravity on the ellipsoid. A positive <5N means that the geoid 
is above the cogeoid. 

The potential of gravity, W is the sum of the potentials of gravi- 
tational force, V and centrifugal force, $ , i.e., W = V + 4> . Now 
we are interested in the change of gravity potential, 6W or the differ- 
ence between the gravity potentials of the masses before and after the 
gravity reduction. Since the potentials of the centrifugal forces of 
the masses before and after the gravity reduction are the same, the change 
of gravity potnetial, 6 W is then the same as the change of gravitational 
potential, 6V . Brun’s formula in equation (2) becomes: 

SN = (3) 

Therefore, the potential change can be referred to as either the change in 
gravity potential or the change in gravitational potential. 

Secondary Indirect Effect on Gravity 

Since Stokes' equation gives the surface of the cogeoid rather than 
the geoid then to be more accurate, the gravity anomaly used in the equa- 
tion must be the boundary value on the computed surface, i.e., the gravity 
station is further brought down from point P Q on the geoid to point 
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Figure 1 : Geoid, cogeoid and ellipsoid 


P on the cogeoid (see Figure 1) before Stokes' equation can be applied. 
In order to bring the gravity station from the geoid to the cogeoid, 
it is done by a simple "free-air reduction" (Heiskanen and Mortiz, 1967, 
p. 142). The condensation reduction is unnecessary because the masses 
between the geoid and cogeoid are relatively small compared to the masses 
above the geoid. The movement of such small mass through the very small 
distance will produce no significant change in SW . Therefore, it can 
be thought that instead cf moving the masses inside the geoid in the 
first place, we immediately put them inside the cogeoid with 6W re- 
maining unchanged (Bomford, 1971, p .’ 493) . 

The change of gravity by free-air reduction in bringing down the 
point for the vertical distance h is computed by 


In practice, the vertical gradient of gravity 3g/9h , is replaced by the 
normal gradient 9y/3h so we have 



= + 0.3086h mgal (4) 

for h in meters. From the above equation, the change in gravity in 
bringing down the gravity station from point P Q on the geoid to point 
P c on the cogeo id can be computed from: 

6 = 0.3086 t$N mgal (5) 

where 6N is the separation between the geoid and cogeoid or the indirect 
effect on geoid undulation in meters. The effect of 6' is called "a 
secondary indirect effect on gravity" (Heiskanen and Moritz, 1967, p. 142.) 

Considering these two indirect effects caused by the gravity reduc- 
tion, the undulation of the geoid, N is 

N = jjj. / f (Ag + 6) S(<p) da + SN (6) 

where Ag is 'the gravity anomaly at point P Q on the geoid. 

Gravity Reduction 

There are several ways in which the gravity reduction can be per- 
formed. Individual methods differ, depending on how the topographic 
masses above the geoid are treated. For example, the Bouguer reduction 
completely removes the topographic masses. The free-air reduction, in 
concept, ignores the existence of the masses so the gravity on the topo- 
graphic surface is reduced to the geoid using the vertical component 
of the gravity gradient. The second method of Helmert's condensation 
shifts and condenses the topographic masses on the layer of the geoid 
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while an isostatic reduction distributes the masses to a certain depth V 

under the geoid. More methods and more details can be found in standard 
textbooks in physical geodesy, e.g. Heiskanen and Moritz (1967), Bomford 
(1971), and Grushinsky (1969). 

Theoretically, all gravity reductions are equivalent and lead to 
the same geoid provided that they are properly applied, i.e., the indirect 
effect is taking into consideration (Heiskanen and Moritz, 1967). In 
practice, however, gravity reductions which give a large indirect effect 
are avoided in determining the geoid. Usually, a free-air anomaly at 
sea level is used in Stokes' equation and is said to be an approximation 
of Helmert's gravity anomaly (Heiskanen and Moritz, 1967, p. 146), The 
indirect effect due to Helmert's reduction is less than 50 cm for the 
topograhic height of 3 km (Heiskanen and Moritz, 1967, p. 145; Grushinsky, 

1969, and Neequaye, 1975); therefore, it may not necessarily be considered, 

However, the required accuracy of geoid undulations may be on the order 
of 10 cm or less in the future. This means that the indirect effects 
cannot be neglected anymore. 

Since the free-air anomaly is mostly used in practice and closely 
related to the anomaly of the second method of Helmert's condensation. 

It is interesting for us to look into more details of these two anomalies. 


Free Air Anomaly vs. Helmert's Anomaly 

Helmert introduced two methods of condensation. In the first method, 
the topography is condensed on a parallel surface located 21 kilometers 
below the geoid (Heiskanen and Moritz, 1967, p. 145), The value of 21 
kilometers, the difference between the semimajor and semiminor axes of the 
earth, is adopted in order to avoid problems connected with the convergence 
of the spherical harmonic series for the potential outside the earth. The 
first method is not popular and not used by geodesists. The second method 
of Helmert's condensation is therefore the one which is referred to in most 
literature and also in this study. 
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In the second method of Helmert's condensation, the topographic f 

y 

masses of volume density p above the geoid are shifted and condensed 
to be a surface layer of surface density ph , where h is the height 
of the topographic surface above the geoid. After condensation, the 
gravity force at point P on the earth's surface is transferred to the 
corresponding point, P Q , on the geoid using the vertical gradient 
of gravity. In a free-air reduction, the gravity force is just 

transferred from point P on the surface to point P Q on the geoid 
by the vertical gradient of gravity, without doing anything with the 
topographic masses. This means that if gravity g obs is measured at 
point P on the topogrpahic surface, the gravity at point P Q on the 
geoid becomes g Qbs + F ‘'here F is the change in gravity by the free- 
air reduction computed by equation (4). After subtracting normal gravity, 
y , at the corresponding point on the ellipsoid, we get a free-air anomaly 
on the geoid as: 

" ®obs * F " V < 7 > 

A gravity anomaly of the geoid obtained by the second method of 
Helmert's fcnndensation, shortly called here as Helmert's gravity anomaly, 

4§ h „ differs from the free-air anomaly at the same point by the amount 
of terrain correction C (Lambert, 1930, p. 116), that is, 

Ag^ * Ag + C (8) 

Equation (0) is obtained based on the linear approximations for the effect 
of topography (Moritz, 1968). We will see later that C is only one 
part of the total attraction change due to the second method of Helmert's 
condensation. 

Considering the procedure in the second method of Helmert's condensation, 
we see that this method maintains the same earth's mass and the same location 
ot the center of the mass of the earth (Bomford, 1971, p. 505). That 
is the second method of Helmert's condensation fulfills the second and 
the third requirements of Stokes' equation. The first requirement is 
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always fulfilled by keeping the definition of the geoid which is the 
equipotential surface having potential . Realizing these facts, 
we see that using the gravity anomalies obtained from the second method 
of Helmert's condensation in equation (6) should give the more correct 
geoid undulations than using the ones from the free-air reduction, By 
using the gravity anomaly of equation (8) in equation (6), the more accurate 
geoid undulation becomes 

N = 7 ^// (Ag + C+6) S(l/>) do + 6N (9) 

= // Ag S(iji) do + jj c S(if>) do 

+ JfS S(iC)do + «N (10) 

N = Nj + ANi + <5N 2 + 6N (11) 

We see that Ni is the undulation which has normally been computed from 
free-air anomalies and used as the 'correct' undulation. This means 
that if the free-air anomalies, Ag are used in Stokes' equation (1), 
there are three kinds of errors incurred in determining the geoid. One 
comes from neglecting the terrain effect in the gravity reduction and 
the other two from the indirect effects. In sequel, we will refer correc- 
tion ANi which comes from the attraction change in Helmert's condensation 
method as an attraction change effect on geoid undual tion; correction 
<5N which comes from the potential change as an "indirect effect" and 
correction 6N 2 as a "secondary indirect effect". 

In the following sections, we will investigate how large each cor- 
rection is. 
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2. Investigation of Formulas for Computation 


It is shown in the previous chapter (equation (11)) that the accurate 
geoid undulation can be computed by: (i) using the gravity anomalies 

obtained from the second method of Helmert's condensation in Stokes' 
equation and (ii) taking the indirect effects on geoid undulation and 
on gravity into consideration. In other words, if free-air anomalies 
instead of Helmert's anomalies are used in Stokes' equation, we must 
compute the corrections caused by the potential change and the attraction 
change inth'e process of Helmert's condensation reduction (the second 
method). In literature, various models of topographic masses are used 
to derive or to estimate the effects of potential change and attraction 
change in the reduction process, for example: a spherical Bouguer plate 

(i.e., a portion of a spherical shell); a plane circular plate 

(i.e., a plane plate with finite radius); and a plane Bouguer plate (i.e., 
an infinite plane plate). These examples of topograpgic models are referred 
to as "the regular part" or "the smooth part" of the topography because 
the irregularity of the topography is not taken into consideration in 
these models. When the irregularity of the topography is considered 
in the model, the masses above and the non-existing masses below the upper 
surface of the regular part of the topography are referred to as the 
irregular part of the topography. 

In this chapter we investigate the existing expressions to see which 
type of model and which expressions we should use to get the geoid compu- 
tation of 10 centimeter accuracy. Before we pick up the appropriate 
expressions for investigation, there is a little confusion to be made 
clear. The confusion concerns the location of the computation point 
where we compute the changes in gravity potential and gravity attraction 
due to the reduction process. The following is a literature review that 
involves the computation point for the potential change and the attraction 
change in the second method of Helmert's condensation. 
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In the derivation of Moritz (1968), he derives the gravitational 
potential of the topographi masses, V , at a point P on the earth's 
surface. The potential of the condensed layers, V s , is implied to 
be lying on point P Q at sea level (Ibid., p. 18). The attraction of 
the topographic masses and the condensed layer, A and A $ respectively 
are both referred to the same point P . By linear approximations for 
the topographic effect, Moritz shows that V = V and A = A -C where 

b b 

C is the quantity known as a terrain correction (Ibid., eq. (67)). 

The change in gravity due to the condensation of topography to sea 
level by Bomford (1971) is evaluated at point P on the earth's surface 
This gravity change is zero if the topography is assumed to be a spherical 
shell or an infinite plane as Bomford states that "Condensing the plateau 
to below the geoid causes no change in g , whether the plateau is regarded 
as an infinite plane or as a spherical shell." There is no definite 
statement, however, concerning the point where the change in potential 
due to the condensation method is evaluated. 

Neequaye (1975) uses a plane circular plate with finite radius 
a thickness b and volume density p as a topographic model (see figure 
2). He obtains a general expression for computing the gravitational 
potential of the topographic masses at some point of height c directly 
above the center of the bottom surface of the circular plate as (Ibid., 
eq. (3-1)): 

V = Trkp[ (c-b) 2 ~ c 2 -(c-b) /a 2 +(c-b) 2 + c /a 2 + c 2 

-a 2 £n(c-b+/a 2 +(c-b) 2 ) + a 2 Jln(c+/a 2 +c T )] (12) 
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Figure 2 : A topographic model used by Neequaye (1975) 


The potential change is the difference between the potential of the actual 
topographic mass and the potential of the condensed layer. To get the po- 
tential of the topographic mass, Vj , the distance c is set equal to b. 
Equation (12) becomes: 


V T = Ttkp (-b 2 +b/a z +b 2 + a 2 £n ) 

I a 


(13) 


This means that the potential of the topography is evaluated at point 
P on the earth's surface. 
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To get the potential of the condensed layer, the quantity p in 
the equation (12) is set to g. By taking the limit b+0, equation (12) becomes 

V = 2-rrkpb (/a 2 + c z -c) . 

V 

Then by setting c equal to b , we get 


V c = 27rkpb (/a 2 + b 2 - b ) (14) 

The potential of the condensed layer is again evaluated at the point 
P . Therefore, the potential change: 

sv - v T - v c 

= Ttkp{b 2 - b/a 2 +~b + a 2 «n ^| . -t ■ l T . } (15) 

is evaluated at the point P on the earth's surface. 

Grushinsky (1969) uses a plane circular plate with an infinite radius 
for computing the potential change. His derivation implies that the 
potential due. to topography refers to point P on the earth's surface 
while the potential due to the condensed layer refers to point P Q on 
the geoid. However, because of the symmetry of a plane circular plate, 
the gravitational potential due to the topography at point P on the 
surface andat point P 0 on the geoid are exactly the same. Therefore, 
the approximate formula for potential change that has been obtained can 
be referred to a point P Q on the geoid. 

Lambert (1930, p. 116) makes a clear statement in regard to finding 
the reduced gravity at point P Q on the geoid. Apart from the free- 
air reduction, we must apply a correction to the measured gravity at 
P on the earth surface. This correction is equal to the attraction 
of the condensed topography at P Q on the geoid minus the attraction 
of the topography at P on the earth's surface. When Lambert derives 
the potentials in isostatic reduction, the compuation points he refers 
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to are the same as those to which he refers for the attractions (Ibid., 
p. 150). Helmert's condensation reduction can be viewed as a special 
case of isostatic reduction where the compensation depth is zero (Heiskanen 
and Moritz, 1967, p. 145 and Bomford, 1971, p. 497). This means that 
the gravitational potential of the topography refers to the point P 
on the earth's surface and the potential of the condensed layer refers 
to the point P Q on the geoid. 

Baeschlin (1948, Chapter XII) uses the point P as the computation 
point for the attraction difference between the actual topography and 
the compensated topography in the isostatic reduction method. He uses 
the point P Q on the geoid as the computation point for potential change. 

Heiskanen and Moritz (1967, Chapter 3) refer the computation point 
to point P on the surface when the attractions of the topography and 
the condensed layer are determined in Helmert's condensation reduction. 

They also point out that to compute the potentials of the topography 
and the condensed layer, the point to which the potentials refer is always 
the point P Q on the geoid. 

From the above literature review, we see some confusion in using 
computation points for the potential change and the attraction change 
in the second method of Helmert's condensation between point P on the 
earth's surface and point P Q on the geoid. To understand what would 
be the reason which causes the confusion, we take a spherical shell to 
represent the topography of the earth. Let the thickness of the shell, 
hp represent the height of the topography and M s be the masses of 
the spherical shell with constant volume density p . The mass M s can 
be computed from: 

M g = p x volume of the spherical shell 

- P*t' it [( R + h p ) 3 " R3 ] 

= If (3R 2 h p + 3Rh 2 + h^) (16) 
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Introduce u = h /R , then 

r 
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M = 47rpR 2 h n (1 + u + fu 2 ) 

bp 3 

When evaluation point at P on the earth's surface. 


To evaluate the potential of the spherical shell at point P on 
the earth's surface, the whole masses of the shell are acting like a 
point mass at the center of the shell because of the symmetry of the 
shell. That is the gravitational potential of the shell at point P 
symbolized by V'ip can be written as: 


V' = — 

|P (R+hp) 

kM s 

R(1 + u) 

Expanding (l+u)" 1 = l-u+u 2 -u 3 + ... and substituting for M s , we have 


V'jp = k*47rpRhp (l+u+-|-u 2 ) (l-u+u 2 -u 3 + ...) 

= 4nkpRh n (l+-~u 2 “4" u3+ •■•) 

U a 3 

To get the vertical attraction of the spherical shell at point P on 
the earth's surface, we have 



< R V 

-kMs 


R 2 ( 1+u) 2 

since (1+u)" 2 = l-2u+3u 2 -4u 3 + ... , then we get 
A'| p = -47rkph p ( l-u+-|-u 2 --f-u 3 + ...) 
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When the spherical shell is condensed to be a layer surface on the 
geoid, this layer has the same mass M s and can still be considered 
as the point mass at the center of the sphere, i.e., the potential of 
the condensed layer, at point P on the earth surface is: 

V ■ ■ v '|p < 23 > 

The vertical attraction of the condensed layer at point P , Agj p is 
also equal to A 'j p . 

When evaluation point at Po on the geoid. 


To evaluate the potential of the spherical shell at point P 0 on the 
geoid, we cannot consider the shell as a point mass anymore since the 
evaluating point is now inside the shell. For a point inside the shell 
at a distance r from the center, the potential of the shell is 
(MacMillan, 1958, p. 38): 

V' = 4?rkp [i(R+h p ) 2 - ~ - -f-r 2 ] (24) 

At point P Q on the geoid r=R , so we have 

V'| p = 4Trkp [*{R+hp) 2 - y r2 " "H 2 ^ 

= 47rkp [±R 2 + Rh p + ih 2 - 1R 2 ] 

= 4rrkp (Rh p + ±h 2 ) 

= 47TkpRh p (lHu) (25) 

Then we have the vertical attraction of the spherical shell at point 
P 0 as: 


l p o = ' ** l p 0 

= ‘ ST (^TrkpRh + 2Trkph 2 ) 


47ikph r 
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For the potential of the condensed layer at point P 0 on the geoid, 
we consider that point P 0 is just outside the condensed layer so the 
layer acts like a point mass, i.e., 



= 4iTkpRh (l+u+-| u 2 ) (27) 

The vertical attraction of the condensed layer at point P 0 on the geoid 
is: 




= - 47rkphp (l+u+|- u 2 ) 


(28) 


A summary of expressions for computing the potentials and the ver- 
tical attractions of the spherical shell and the corresponding condensed 
layer at point P on the earth's surface and at point P 0 on the geoid 
is given in Table 1. Also shown in the table are the changes in potential 
and vertical attraction at points P and P 0 . 


Comparing the quantities of V* , V' , A' and A' at point P 
on the earth's surface and point P Q on the geoid, we see that the dif- 
ferences are in the order of u . Neglecting u-terms in the computation 
would cause a relative error smaller than 0.14% (Moritz, 1968, p. 4). 

To get an idea how -large the difference of these quantities can be, we 
take an example of hp = 3000 meters, R = 6371000 meters, Newtonian gravita- 
tional constant k = 6.672 x 10 -11 m 3 /(kg-sec 2 ) and volume density 
p = 2.67 gm/cm 3 or 2.67 x 10 3 kg/m 3 . The result is shown in Table 2. 

From Table 2, using different evaluation points between P and P 0 
gives differences of potentials in the order of 1-2 kgal-m and differ- 
ences of vertical attractions in the order of 0.3-0. 6 mgal for hp = 

3000 m . By Brun's formula of equation (2) , the potential change of 
1 kgal-m gives approximately the indirect effect on geoid undulation 
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Table 1 ; Potentials and vertical attractions of the masses 
of a spherical shell and the condensed mass onto 
the geoid at two evaluation points. 


Quantities 

Evaluated at P 0 

Evaluated at P 

V 

4irkpl!h p (l+f) 

RTRjj * 4 tt kpRhp( I+- 5 - - y + ... ) 

V s 

4irkp Rhp ( l+u+-j- 2 ) 

47TkpRhp(l+y “ -jj- + ... ) 

V . 

-4TrkpRh p (y+-j ) = 
[ -^kPh p 2 ] 

0 

A' 

-4tt kphp 

(R+h‘% " - 4 TTkPh p (l-u +-3 3 - + *’• > 

r 


-4TTkphp(l+U+y ) 

-4TTkPh p (l-u+^j- - 3 - + ••• ) 

A' - Aj 

4TTkphp(U+y ) 

0 


Note 



= effect due to spherical 
shell 

= effect due to the condensed 
layer 

= Newton's gravitational const. 
= mass of the spherical shell 
= density of the earth's crust 
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Table 2: Numerical example of Table 1 with h n = 3000 m , Y 

R = 6,371 km , k = 6.67 2x 10" 11 m 3 /(kg-sec 2 ) 
and p = 2.67 x 10 3 kg/m 3 . 


Quantity (in unit) 

Evaluated 
at P 0 

Evaluated 
at P 

1 

Difference be- 
tween two evalu- 
ation points 

V' (kgal-m) 

4279.6501 

4278.6430 

1.00’, 1 

Vg (kgal-m) 

4280.6578 

4278.6430 

2.0071 

V' - Vg (kgal-m) 

-1.0077 

0 

-1.0077 

A* (mgal) 

-671.5810 

-671.2650 

-0.3160 

Ag (mgal ) 

-671.8973 

-671.2650 

-0.6323 

A' - Ag (mgal ) 

0.3163 

0 

0.3163 


of 1 meter. Consider that the average terrain is lower than 3000m, we would 
normally have the differences less than the values shown in Table 2. In 
the past, the magnitudes of differences mentioned above can be tolerated. 
Therefore, this could be a reason that not much attention was given to 
seriously define the location of evaluation points for computing the 
potential change and the attraction change due to the second method of 
Helmert's condensation. Now we are talking about the computation of geoid 
undulation in the order of 10 centimeters, therefore, we need to define 
the exact location of the evaluation point. To do this let us follow the 
reduction process in the second method of Helmert's condensation step-by- 
step as follows, see also figure 1: 

(a) We have gravity value, g Qbs at point P on the earth's surface. 
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(b) Remove the topographic masses above the geoid. We compute the at- 
traction, Af at the point P caused by these masses and subtract 
from the observed value, g 0 jj S , 

(c) Bring the masses back and condense on the layer of the geoid. Here, 
we have to compute the attraction, A c at the point P again which 
is caused by the condensed layer and add to the result in (b) . 

(d) Transfer the gravity to point P 0 on the geoid by a free-air reduc- 
tion, F . For the point P above the geoid, F is added to the 
result in (c) . 

(e) Compute a normal gravity, y of the reference ellipsoid at corresponding 
point Q 0 on the ellipsoid. 

(f) We finally get 

Ag H = 9 obs " A T | P + A C|P + F “ Y |Q 0 ^ 

The vertical line followed by a letter indicates the point where 
a particular function is evaluated. 

The secondary indirect effect is assumed to be negligibly small 
and not included in the reduction procedure explained above. We can 
see from equation (29) that g^* F - y is in fact a free-air anomaly 
at the geoid and (Ay- Ac ) is the gravity attraction change due to the 
second method of Helmert’s condensation and it is evaluated at point 
P on the earth's surface. 

To determine the potential change in the reduction process, we realize 
that the indirect effect, <SN on the geoid undulation is caused by the 
potential change, 6W on the geoid due to the change of the position 
of the masses in the reduction process. Then we have to evaluate the 
potential change, SW at point P 0 on the geoid. 

We can conclude that the computation point or the evaluation point 
of the attraction change is at point P on the earth's surface and the 
computation point of the potential change is at point P 0 on the geoid. 
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In the next section we wfl'< explore the change in gravitational 
potential due to the second method of Helmert's condensation, which will 
lead us to the indirect effect on geoid undulation. 


2.1 Effect of the Second Method of Helmert's Condensation on Potential Change 

The gravitational potential of the topographic masses at point P 0 
is: 


v = k ///t 5 - "" (30) 

The integral is extended over the masses outside the geoid; dm is the 
element of mass and l is the distance of mass dm from the evaluation 
point P Q , (See Figure 3). Let p be a volume density of the topo- 
graphic masses and assumed to be a constant. Then equation (30) can 
be written as: 

V = kp/JJ (31) 

where dv is the element of volume and the integral is extended over 
the exterior <?f the geoid. The element dv = r 2 dgdr where da is the 
elements of solid angle and r is the distance between mass element 
dm and the center of the earth. By planar approximation, 

dv = R 2 dadz 

(see Figure 3) so that equation (31) becomes 

v . kpR 2 // / h 

a z-0 K 

where a denotes solid angle over the whole sphere. Wc rewrite this 
equation again as 
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V . kpR 2 II f hp ♦ kpR 2 11 I h ^ 
o 2*0 0 Z*hp 

= v' -f V" 


Oilimu V:^i2 : , 

OF POOH QUALITY 


( 32 ) 


We see that V' represents the gravitational potential at point P 0 
of the plate which is called the regular part of topography. V" repre- 
sents the gravitational potential at point P Q of the masses above the 
plate and non-existing masses under the plate which are called the ir- 
regular part of topography or topographic variation. That is the potential 
of the topographic masses can be split into the potentials of the regular 
and irregular parts. 


Similarly, the potential of the condensed layer after Helmert's 
condensation can also be split into two corresponding parts: 

V s = V,; + VJ (33) 


Q 



The regular part of topography is represented by 
a spherical Bouguer plate and the irregular part is 
represented by the masses above the plate and the 
non-existing masses under the plate (cross hatched 
areas in the figure). 

- 20 - 


ORIGINAL PAGE IS 
01- POOR QUALITY 


9 

Y 


Then the potential change at point P 0 , <$V in Helmert's condensation 

reduction can be written as 

6V = V - V c (34) 

= (V 1 + V") - (V£ + Vp 
= (V' - V*) + (V" - vp 

= 6V + 6V" (35) 

where <$V' represents the potential change of a regular part of topo- 
graphy and SV" represents the potential change of an irregular one. 
Therefore, we can investigate the expressions for the potential changes 
of a regular part and an irregular part of topography separately. 

2.1.1 Regular Part 

As mentioned before, the effect of a regular part of topography 
or the Bouguer plate can be computed by different degrees of approximation. 
We start from the one that has the least degree of approximation. 

(a) Spherical Bouguer Plate : The regular part of the topography 

is described as a portion of the spherical shell with thickness hp and 
angular radius \p from the computation point. The gravitational potential 
of this topographic mass at point P Q on the geoid can be written as 
follows (Baeschlin, 1948, p. 519) and (Lambert and Darling, 1936, p. 6): 

(36) 


( 37 ) 
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V s ukpR 2 {sin 2 ijj [f(i|>, h p ) - f(ip, 0)3 + 4B Q } 

where 

f(*,x) = -SSiiUt- + cos » tan Y + C os* in tan (4- + i) 
3cos 3 y cosy 4 c 

with 

tany= tanf + 


) . .. •*4**^| 
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The angle ip is the angular radius of the plate subtended at the center 
of the sphere (see Figure 4), and B Q is defined below. 


The gravitational potential of a spherical Bouguer plate can also 
be expressed in hyperbolic functions as (Baeschlin, 1948, eq. 75.27; 
and Lambert and Darling, 1936, eq. 28c): 


V = 4'nkpR 2 {B 0 +-|-sin 3 'Mcosh 3 <i>2 - cosh 3 <h) 

+ ■7 _ sin 2 ^cos , M4>2-<f>i) 

+ -^-sin 2 ^cos^( sinh 4>2 coshcf> 2 ~ s1nh$i cosfvh)} 

with 

B„ * " *4) Cl + f (>)] 

sinh 4 >x = tan 

h 

sinh * 2 - tan f + * 5 ^ 



Figure 4 : Topography is represented by a spherical Bouguer plate 


To obtain the gravitational potential of the condensed layer, we 
make use of the expression for the potential of a spherical disk given 
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by Baeschlin (1948, equation 75.17) as: 

V' = 2rrkpdr (E-z) 

with E = /R * + r 7 - 2Rrcosi|J where r is the distance from the center 
of the earth to the evaluation point, dr is the thickness of the disk 
and z = r-R . To apply the above equation to the condensed layer of 
a spherical Bouguer plate of height h p , we replace dr by h p and 
with the evaluation point at P 0 on the geoid, we have r=R and z=0 
then, 


E =/R 2 + R 2 - 2R 2 cos ip 
= R /2 - 2 cos £ 


With the trigonometric identity (1-cos<|0 = 2sin 2 we get 
E = 2R sin \ 

So we can write the potential of the condensed layer of the spherical 
Bouguer plate at point P 0 as 

Vg = 47rkph p R sin|- (40) 

Therefore we can compute the potential change at point P 0 of spherical 
Bouguer plate with equation (36) or (38) and equation (40). 

Next we turn to a less accurate way to estimate the potential of 
the regular topography. 


hr 

5 

V 


Figure 5 : Topography is represented by a plane circular plate 
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(b) Plane Circular Plate : Here, the regular part of the topography 
is a flat plate with thickness h p and linear radius a from the compu- 
tation point. The expression for the graviational potential of the plate 
Is less complicated than the pontential of the spherical Bouguer plate. 

The potential of the plane circular plate at point P 0 is given by 
(Heiskanen and Moritz, 1967, p. 128; Baeschlin, 1948, p. 500): 

V' = Trkp[-hp + hpEi + a 2 £n ] (41) 

with ti = /a 2 + h p (see Figure 5). For an equivalent area with the 
spherical plate in (a), the linear radius of the plane plate is: 

a = 2R sin|r (42) 

Using equation (3-9) of Heiskanen and Moritz (1967) and letting 
quantity c which is the height of the evaluation point equal to zero; 
we get the potential of the condensed layer of the plane circular plate 
at point P 0 as: 

Vg = 2irkKa 

= 2?rkphp(2R sin-Sr) 

= 4-nrkph R sin|- (43) 

Comparing equation (43) with equation (40), we see that the potentials 
of the condensed layers of the equivalent size of plates are the same 
regardless of a spherical or a flat plate. 

Finally, we turn to the least accurate approximation for the regular 
topography. 

(c) Infinite Plane Plate or Plane Bouguer Plate : This is the plane 
circular plate of height h p in (b) with a radius of infinity. Grushinsky 
(1969, p. 210) gives for the potential of the infinite plate at point P 0 : 

r = 2rtph p a [1 - Ijf. + 
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and the corresponding potential for the condensed topography is; 


= 2Tikphp • a 


Of specific interest, here is the difference V' - (see eq. 

35) which can be simply computed for model (c) but not models (a) and 
(b). We have: 

SV - -2irkph p • a [£§- - (Jg-} 2 ] 

- - trkphp ♦ iSgL 

As the radius a approaches infinity, the second term in the equation 
goes to zero, so we have 

<5V' = “Trkphp (44) 

Equation (44) will be referred to later as Grushinsky's formula. 

2.1.2 Irregular Part 

The gravitational potential at point P 0 of the irregular part 
of topography is defined to be (equation (32)): 

V" = kpR 2 // / (45) 

a z=h p 

To solve equation (45), we consider the planar approximation of the topo- 
graphy shown in Figure 6. From the figure we have 

Z 2 = Zl + z 2 (46) 


or we can write 


1 . = 1 (i + 4!)-* 
z z 0 v I V 
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Expand the expression in parenthesis as a binomial series, we get 


1-1/1 L Z 2 . 3 Z 4 15 Z 6 


i_-2.fl 1 + 3 £ _1_5 

* " *• u ' 8 *; ‘ * 8 *<j 

6 


1 

To 


*o 


l t- .3 2T 15 r . v 

* - a 8 ^"5"" ~ L n "« 7 * * * / 


+ . . . ) 


48 


Substitute l/£ into equation (45) and integrate to obtain 




Vi + V2 + V3 + V4 + . . . 


(47) 


with 


V. ■ kpR 2 / J / h 


0 z=ht 


= kpR 2 // 


£0 


z=h r 


da 


- kpR 2 ffkte- da 


(48) 


V 2 = kpR 2 // / - i - dz da 


a z=h r 


= - i kpR 2 ff-f- 




h 

z=h f 


da 


h 3 -hr 


= - -g-kpR 2 da 


(49) 


h , z 4 


V 3 = kpR 2 // / -jp — dz da 
a z=h n 


= O k P R2 II 


h 5 -ht/ 

a ^0 


da 


V 4 = 


kpR 2 // j" 
a z=h. 


H- — dz da 


48 £ 0 7 


KpR 2 II . 

a *0 


h 7 -hp 7 


da 


(50) 


(51) 


etc . 
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Figure 6 : Topography in the planar approximation. The geoid 

may be visualized as a flat plate (Moritz, 1966, p. 33). 
The regular part of the topography is then visualed as 
the plane Bouguer plate. 


Equation (47) with equations (48) to (51) resembles equation (22) 
of Moritz (1968) where he derives the gravitational potential at point 
P on the earth's surface of the irregular part of topography: 

V"| p = V* + V* + V* + . . . (52) 

where 

Vf ■ kpR 2 // da 
a 

V* = - '4-k PR 2 // da (53) 

O 



etc . 


The distance can be seen in Figure 3 on page 20. 

The potential of a surface layer at point P 0 on the geoid is given 
by (Heiskanen and Moritz, 1967, eq. 1-16): 

V s = R 2 da (54) 
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where k is the surface density of the layer. Applying this equation 
to the irregular part of the topography, we have surface density k = 
p(h-hp) . Then the potential of the condensed layer of the irregular 
part at point P Q is; 

V" = k ;j R* do 

s a u 

- kpR z // da (55) 

a lo 

We see that VJ is equal to Vi term of equation (47) so that <$V" 

(=V"-V^) is: 

6V" = V 2 + V 3 + V 4 + . . . (56) 

Therefore, equation (56) in conjunction with the potential change 
of either model (a), (b) or (c) in subsection 2.1.1 can be used to obtain 
the total potential change due to the second method of Helmert's conden- 
sation. 

In the next section we will look at the effect of the second method 
of Helmert's condensation on the attraction change. 


2.2 Effect of the Second Method of Helmert's Condensation on Attraction 
Change 


By definition, the vertical gravitational attraction at point P 
on the surface is written as: 


= . 

3r p 


where V| p is the gravitational potential at point P and r p repre- 
sents the vertical direction at point P for the spherical approximation 
of the geoid. Comparable to equation (32), we may separate the gravi- 
tational attraction of actual topography into two parts. Since V = 

V' + V" , then 
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a=A' + A" 


f\? s r? |» 

POOR QUAINT 


9r r 


Similarly, the attraction of the condensed layer at point P can be 
written as 


A s = A* + Ag (58) 

Then the difference of gravitational attraction between the topography 
and the condensed layer at point P is: 

Ay - Aq “ f)(\ = A “ Ag 

= (A' + A") - (A' + a;) 

= (A' - A^) + (A" - AJ) 

= 6A* + <5A" 

Therefore <sA 1 represents the attraction change of a regular part of 
topography and fiA" represents that of an irregular one. So we will 
investigate the attraction change of the regular part of topography sep- 
arate from the attraction change of the irregular one. 


(59) 

(60) 


2.2.1 Regular Part 

The attraction change of the regular part may be computed in three 
different ways as follows: 

(a) Spherical Bouguer Plate : The attraction of the topography 

at point P on the earth's surface given by Baeschlin (1948, p. 529) is 


A' = -2rrkpR s i nip{ (cosh4> + sinh0) (■£• costp + sin 2 ^ s i nh 2 ^) 
- -§- sin 2 ijj coshtj) - simfj cosiptj) }^ 2 


(61) 


with 


sinh 


r-RcosJ; 

Rsini 


4>i : r = R 

<f >2 : r = R+hp 


( 62 ) 
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The expressions for si nh^i and sinh <#>2 in equation (62) can be easily 
reduced to the corresponding expressions in equation (39). The attraction 
of the corresponding condensed topography at the point P is 

A s * I 1 - ) (63) 

where 

E = /TR+h p ) z + W - 2R(R+h p ) cosijJ (64) 

The attraction change of a regular part of the topography is then calculated 
from 


9 

Y 


<5A' = A' - A' . 

j 

(b) Plane Circular Plate : The attraction of the actual topography 
at point P on the earth's surface given by Heiskanen and Moritz (1967* 
eq. 3-6) is: 


A' = 27rkp (a + hp - /a 2 + h p ) 


(65) 


and the attraction of the condensed topography at the point P is (Ibid., 
eq. 3-10) : 


A' = 27rkph n (1 


/a 2 + h p 


) • 


( 66 ) 


(c) Infinite Plane Plate : The attraction change of the infinite 
plane plate or Bouguer plate is zero. This can be verified by equations 
(65) and (66) taking a . We can write equation (65) as 


A' = 27rkp (h p + a - a A + (^-) 2 ) 


With a oo , -f- + 0 and 


A' = 27rkph r 


/a 2 + hp 


■* 0 , the above equation becomes 


and equation (66) becomes 
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Ag “ 2irkph p . 
Therefore, 

<5A' = A* - Ag B 0 . 
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Next we turn to the attraction change of the irregular part of topo- 
graphy. 


2.2.2 Irregular Part 

The total attraction change (sum of the attraction change of the 
regular part and that of the irregular part) at point P on the earth's 
surface in the second method of Helmert's condensation is given by Moritz 
(1968, eq. 63) to be: 

6A = A - A = - C (67) 

where C is the quantity called the "terrain correction". By the planar 
approximation (Ibid., eq. 42), 

C = i kpR 2 // ' h ~ h P )2 do (68) 

o H 0 3 

In his derivation, Moritz uses a spherical shell as a regular part of 
topography (Ibid., p. 6). If we examine the derivation, we can see that 
the total amount of attraction change (-C) comes from the attraction 
change of the irregular part alone. This argument is strengthened by 
the fact that the attraction change of a spherical shell at point P 
on the earth's surface is zero (see Table 1 on page 16) i.e., we have 

<5A" = - C (69) 

It should be noted here that the masses further away from the eval* 
uation point P contribute a very small effect on the terrain correction. 


- 31 - 


ft 


The effect decreases very rapidly as we can see from the term A© 3 
in the denominator of equation (68). 

Using equation (69) and the attraction change of the regular part 
(model (a) or (b) In subsection 2,2.1), we arrive at the total attraction 
change within a certain cap size. This statement does not contradict 
the statement we just made that the attraction change of the regular 
part of topography for the whole spherical shell is zero. In practice, 
we do not integrate over the whole sphere. We usually consider the effect 
(e.g., the effect of gravity anomalies on the geoid undulations) within 
a certain cap size, for example, 10° or 20° from the evaluation point. 

To compute the attraction change we may also have to consider the effect 
of the topography only within the cap size. 

In the next section, we will compute the numerical values of the 
potentials and attractions presented here. Then we will make a comparison 
so that the conclusion concerning which expressions vie should use to 
obtain the most accurate geoid computation can be made. 

2.3 Comparison of Formulas Using Simple Test Mode ls 

As discussed earlier, the regular part of the topography may be 
represented by at least three models, i.e. the spherical Bouguer plate, 
the plane circular plate and the infinite plane plate. We, then., want 
to know whether it makes any difference to use one model instead of another 
and whether there should beliniitations to do so or not. Moreover, based 
on the expressions for computing the potential change and the attraction 
change of the irregular part of the topography, we would like to know 
if these expressions have any relations with the expressions for the 
regular part. 

The topographic masses above the sea level occur in about 30% of 
the total area of the surface of the earth. Masses further away from 
the compuation point have very small effects on the values of the potential 
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change and the attraction change. It would then be more appropriate 
if we used the regular part of the topography at a certain cap size to 
cover a particular continental area only. This is because there are 
no topographic masses above the geoid to be condensed in the ocean area. 

Thus the effects on the potential change and the attraction change at 
any evaluation point computed from ocean area masses are zero. By selecting 
one cap size for the whole area we may however take too much mass into 
account in the computation of the regular part. If it is possible, then 
we would like to see that the effect computed from the irregular part 
will compensate for the excessive effect from the regular part. 

To investigate the numerical effect on the potential change and 
the attraction change of the regular part of topography, spherical Bouguer 
plates with height h varying from 500m to 9000m and angular radii & 
varying from 1° to 180° are used in this study for the expressions of model 

a 

(a) in subsections 2.1.1 and 2.2.1. Plane circular plates with the same 
heights and radii as those of spherical Bouguer plates are used for the 
expressions of model (b) in subsections 2.1.1 and 2.2.1. The same size 
of the plate in each model is used so that we can compare the effects on 
potential change and attraction change. 

To compare the effects oh potential change and attraction change 
of the regular part with those of the irregular part, we use a "spherical 
ring" and a "cylindrical ring" (see Figure 7). A spherical ring is a 
portion of a spherical Bouguer plate which is bounded by angular radii ipi 
and \p 2 while a cylindrical ring is a portion of a plane circular plate 
bounded by linear radii a* and a 2 . 

The potential change at point P 0 of a spherical ring can be computed 
from the difference between the potential changes of spherical Bouguer 
plates of angular radii ipi and t explained in subsection 2.1.1(a). 
Similarly, the attraction change at point P of a spherical ring is 
computed from the difference between the attraction changes of spherical 
Bouguer plates of radii ip i and i p z explained in subsection 2.2.1(a). 

This method of computing the potential change and attraction change from 
two different radius of spherical Bouguer plates will be referred to 
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(a) A section of a spherical (b) A cylindrical ring 

ring 


Figure 7 ; Models for comparison between the effects on the potential 
change and the attraction change of the regular part and 
the irregular part of the topography. 
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later as TDSBP (Two Different Spherical Bouguer Plate) method. 

The potential change at point P 0 of a cylindrical ring can be 
computed from the difference between the potential changes of plane circular 
plates of linear radii ai and az explained in subsection 2.1.1(b). 

The attraction change at point P of a cylindrical ring is computed 
from the difference between the attraction changes of plane circular 
plates of radii ai and az explained in subsection 2.2.1(b). The 
method of computation here will be referred to as TDPCP (Two Different 
Plane Circular Plate) method. 


The potential change of the irregular part of topography (equation 
(56) with equations (49) to (51) in subsection 2.1.2) can also be used 
to compute the potential change of a cylindrical ring as follows: From 
equation (49), 


Vz = - f kp // h -&- R 2 do 

a Hi 


For a cylindrical ring, h is constant and we have R 2 da = rdadr and 
H 0 = r (see Figure 7 (b ) ) . So we can write 


V 2 =--l-kp (h 3 - h 3 ) / az J 2lr da 
v r=3i a=0 

» - (h 3 « hp) J * 2 it 


rdr 

r 3 


■ - i nkp (h 3 - h 3 ) <-±) 


3z 

r=ai 


■ - f (h 3 - h 3 ) (70) 

Similary, V 3 and V 4 terms for a cylindrical ring with height h and 
radii ai , a 2 can be computed directly from 


V 3 = *7 0 Trkp (h 5 - h®) (—-I 1 ) 




a 3 a 3 
a i a 2 


and 


V 4 


O’ 7 - -hj) <4jfr> 


(71) 
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The attraction change of the irregular part of topography is equal 
to the negative value of the terrain correction. Equation (69) with equa- 
tion (68) can also be used to compute the attraction change of a cylin- 
drical ring. By the same algebraic manipulation used for V 2 , V 3 and 
V4 terms, we can write 

SA" - -Ttkp (h-h ) 2 (73) 

The method of computing the potential change and attraction change 
of a cylindrical ring by using the expressions for the irregular part 
of topography explained above will be referred to as an integration method. 

A spherical ring and a cylindrical ring are called identical if 
the linear radius of the cylindrical ring is related to the angular radius 
of the spherical ring by equation (42) and the heights of both rings 
are equal. By using the identical ring in TDSBP, TDPCP and the integra- 
tion method, we can compare the values of the potential changes and the 
attraction changes. Then we may find a relation between the expressions 
of the regular part and the expressions of the irregular part for the 
potential change and the attraction change. 

The next .subsection is the numerical results in computing the potential 
change using the various expressions explaiend above. 

2.3.1 Comparison of Formulas for Potential Change 

The potential change in the second method of Helmert's condensation 
is the difference between the potential of the topography and the potential 
of the condensed layer, evaluated at point P 0 on the geoid. The poten- 
tial change at point P 0 on the geoid of a plane circular plate as computed 
from equations (41) and (43) are negative. Its values become more negative 
when the thickness of the plate increases. It is approximately a parabolic 
function of the thickness or the height of the plate. Figure 8 shows 
the typical relation between the potential change, 6V' and the height, 
hp of a plane circular plate when the radius of the plate is held fixed. 




9 

y 

The numerical values of SV' in Figure 8 are obtained for the radius of 
the plate being 5° and the earth radius, R being 6,371 kilometers. 

The potential change at the point P 0 of a plane circular plate 
is also more negative when the radius of the plate increases and the 
height of the plate is held fixed. However, the value of the potential 
change is changing very slowly in this case. As we already discussed 
in subsection 2.1.1, an infinite plane plate is a special case of a plane 
circular plate where the radius of a circular plate is extended to infinity. 
Subsequently, we should expect that as the radius of a plane circular 
plate becomes larger, the value of the potential change of the plate 
approaches the value computed from Grushins, ky 's formula (equation (44)) 
which is the potential change of the infinite plane plate. This intuition 
is true as we can see from Table 3. Table 3 shows the potential change 
of a plane circular plate in relation to height hp from 500 m to 9000 m 
and to angular radius from 1° to 20°, along with the value computed by 
Grushinsky's formula. Table 4 gives the potential change error in percents 
if we were to use Grushinsky's formula to compute potential change instead 
of using a certain size of a plane circular plate. For example, if we 
use Grushinsky's formula for the plate of height 3000 m, the potential 
change is -0.504 kgal-meter or the indirect effect on geoid undulation 
is about -50 cm. Suppose that the plane circular plate of 5° radius 
should be used instead of using Grushinsky's formula in the computation. 

We have the relative error of 0.2% which is equal to 0.1 cm in computation. 

In this particular example we see that the difference between the two 
expressions is insignificant for the geoid accuracy of 10 cm. Even in 
the extreme case when the height of the terrain is 9000 m above the geoid, 
the error in using Grushinsky's formula to compute the potential change 
at point P 0 is less than 3 cm. We can, then, say that Grushinsky's 
formula may be used as an approximation to compute the potential change 
at point P 0 of the plane circular plate. 

The potential change of a spherical Bouguer plate at the point P Q 
is compute from equations (36) and (40). The value of the potential 
change of a spherical Bouguer plate starts with nearly the same negative 
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Table 3 : (continued) 
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Jable_4: Potential change error of Grushinsky's formula (equation (44)) 

when compared to result for a plane circular Dlatp fin riovron 
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Table 4 : (continued) 
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value as that of a plane circular plate for the radius of 1°. The po- 
tential change of a spherical Bouguer plate however is less negative 
with the increase of radius ip . The value becomes zero around = 

39° , and beyond this point, the potential change is positive. As an 
example, Figure 9 shows the potential changes of the plate 3OQ0 m high 
with respect to radius ip computed from equations (36) and (40) for 
a spherical Bouguer plate, equations (41) and (43) for a plane circular 
plate and equation (44) for an infinite plane plate. Also shown in the 
Figure 9 is the difference of potential changes between the spherical 
Bouguer plate and the plane circular plate. The difference shows almost 
a linear function of radius t p . 

Comparing the potential of a spherical Bouguer plate computed from 
trigonometric functions (equation (36)) and that from hyperbolic functions 
by equation (38), it is found that both equations give the same numerical 
value up to 10- 6 kgal-meter. 

The computation of the potential change of the irregular part of 
topography, 6V" by equation (56) with equations (49) to (51) shows 
that V4 term has a very small contribution to 6V" . For example, 
using the cylindrical ring explained in Section 2.3 with radii between 

* 

1° to 2° and height 6000 m, the contribution of V4 term is just 10" 6 
of <$V" . The contribution of V4 relative to SV" is even smaller 
when the radii of the ring are large and/or the height is lower. There- 
fore, in most cases we do not need to use Vi f term in the computation 
of <5V" . 


Comparing Two Different Spherical Bouguer Plate (TDSBP) method, 
Two Different Plane Circular Plate (TDPCP) method and the integration 
method for determining the potential change at point P 0 on the geoid 
of the ring of width 1° (the difference between angular radii and 
\pz or their equivalent linear radii a x and a 2 ), a good agreement 
was found between TDPCP method and the integration method. This means 
that the expressions for potential change of a plane circular plate 
(equations (41) and (43)) is compatible with that of irregular part of 
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POTENTIAL CHANGE 
(kgal -meters) 


9 


okknal h 
of poor qimuyv 



Figure 9 : Potential changes of the regular part from different 

expressions for a spherical Bouguer plate, a plane cir- 
cular plate and an infinite plane plate. 
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9 

topography (equation (56)). Table 5 is the example showing the comparison ¥ 

TDPCP and integration method for the height of 3000 m. The discrepancy 
is less than 10“ 9 kgal -meters which is negligile. 


2.3.2 Comparison of Formulas for Attraction Change 

The attraction change in the second method of Helmert's condensation 
is the difference between the vertical attractions of the topography 
before the condensation and the condensed layer, evaluated at point P 
on the earth's surface. The attraction change is subtracted from the 
free-air anomaly to obtain the Helmert's anomaly. Table 6 and Table 7 
show the attraction changes at point P on the surface when the topography 
is represented by the spherical Bouguer plate (SBP) using equations (61) 
and (63) and the plane circular plate (PCP) using equations (65) and 
(66). Height hp of the plates varies from 500 m to 9000 m and spherical 
radius, 1 1 > , varies from 1° to 20° . Comparing the values of the at- 
traction changes of SBP and PCP for small plate (say < 10°) reveals 
that the attraction change of SBP is about three times larger than that 
of PCP. 

The attraction change of PCP becomes smaller as the radius of the 
plate increases and vanishes when the radius is infinitely large as we 
have proved in subsection 2.2.1(c). The attraction change of SBP com- 
puted from equations (61) and (63) however does not go to zero when 
^ = 180° or when the plate becomes the whole spherical shell. The ex» 
pected value of a spherical shell at the point P is expected to be 
zero (see Table 1 on p. 16). Figure 10 shows the behavior of the attrac- 
tion change of SBP with height 3000 m as the spherical angle increases 
from 10° to 180°. The minimum value of the attraction change occurs 
approximatley at ip - 70°. This minimum point is compatible with the 
curve of the potential change of SBP plotted in Figure 9 where the maximum 
slope of the curve is also approximately at ^ = 70°. This is true regardless 
of the evaluation points for the potential change and the attraction change 
not being the same. 
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Table 5 : Comparison of the potential change of a cylindrical ring 
3,000 m high between Two Different Plane circular Plate 
(TDPCP) method and integration method, ih is the inner 
radius of ring and ai is its corresponding linear dis- 
tance taking R « 6371000 meter. 6V' is the potential 
change computed by the expressions for the plane circular 
plate and 6V" is that for the irregular part of topography. 


1 Pi 

(deg.) 

ai 

(m) 

<5V" 

(kgal-m) 

<5V 

(kgal-m) 

<$V" - 6V 
(kgal-m) 

i 

111194. 

1 

0.226400700-02 

0.226400700-02 

-0. 1 10600400-12 

2 

222879 . 

0.754050180-00 

0.784050150-00 

-0. 102521 120-11 

o 

000547 . 

0.077411400-00 

0.07711 1400-00 

0. 120406 180-12 

4 

4448119 , 

0,226426480-00 

0. 226426450-00 

0.801091720-10 

8 

885790. 

0. 150902620-00 

0, 150902620-00 

0. 121924060-10 

O 

660085. 

0, 107792900-00 

0. 107792910-00 

-0,802166700-11 

7 

777UU0. 

0.000004250-04 

0.000004210-04 

0,079978910-11 

O 

UU0007 . 

0.620550710-04 

0.620550700-04 

-0.647800400-11 

9 

999726. 

0.502727940-04 

0.502720000-04 

-0,080174810-11 

10 

11 10580. 

0.411210050-04 

0.41 1217000-04 

0. 197802000-10 

u 

1221267. 

0.042505400-04 

0.042505500-04 

-0. 180841090- 10 


1001902. 

0.209790960-04 

0 . 209790900-04 

0.61 1800990-1 1 

io 

1412405. 

0.240009450-04 

0.240009710-04 

-0.285066600-10 

14 

1552059. 

0.215 124110-04 

0.215120900-04 

0.210476480-10 

18 

1660166. 

0. 100 160900-04 

0. 10!! 16 1050-04 

-0. 145292160- ! 0 

16 

1770044, : 

0. 165985700-04 

0. 165955700-04 

0.800802690-11 

17 

1000000. 

0. 147451400-04 

0. 147481400-04 

-0. 1 10477600- 12 

10 

1990200. 

0. 101000640-04 

0. 101060720-04 

-0.794092660-1 1 

19 

2100007. 

0. 110620100-04 

0. 110620040-04 

-0. 162690100-10 
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Table_6: Attraction change due to condensation when the topography is represented 

by a spherical Bouguer plate (in mgal) as computed from equations (61) and (63). 
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Attraction change due to condensation when the topography is represented 
by a plane circular plate (in mgal) as computed from equations (65) and (66) 
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Table 8 ; Comparison of the attraction change of a cylindrical ring 

3,000 m high between Two Different Plane Circular Plate (TDPCP) 
method and integration method. t|»i is the inner radius of ring 
and ai is its corresponding linear distance taking R = 
6371000 meter. 6A ' is the attraction change computed by the 
expressions for the plane circular plate and 6A" is that 
for the irreglar part of topography. 


*1 

(deg.) 

a 1 
(m) 

SA" 
(mgal ) 

6 A' 
(mgal) 

6A" - <5A 1 
(mgal) 

l 

111194. 

.... . . r . r J 

0 . 2264B290D-02 

0 . 22626576D-02 

0 . 2 1623922D-05 


222379. 

0. 75490 166D-03 

0.754684 16D-03 

0.21749959D-G6 


333547. 

0.37742207D-O3 

0.377369 11 D-03 

0.52953204D-07 


444609 . 

0 . 22643022D-03 

0.22641 136D-03 

0. 18859667D-07 

Bii 

555793. 

0. 15O93429D-03 

O. 15092595D-03 

0 . 83374767D-08 

H 

666065. 

0. 1O779375D-03 

0. 10770951 D-03 

0.4241 145BD-08 

B 

777000. 

0. 0003090 ID-04 

0 . 80828520D-04 

0.2381 3839D-08 

mm 

800037. 

0 . 62855659D-04 

0 . 6285422OD-04 

0, 14309949D-08 

B 

999726. 

0. 50272979D-04 

0 . 50272058D-04 

0. 9203309 8D-09 

10 

1 1 10538. 

0.41 121928D-04 

0.41 121312D-04 

0.6 1583303D-09 

11 

1221267. 

0. 34258629 D-04 

0.3425820 ID- 04 

0 . 42754606D-09 

12 

1331902. 

0. 28979 157D-04 

0.2897885 ID-04 

0 . 30608295D-O9 

13 

1442435. 

0 . 24830990D-04 

0 . 24830765D-04 

0.22491537D-09 

14 

1552859. 

0. 2 15 12445D-04 

0.215 12276D-04 

0. 168925 10D-09 

15 

1663165. 

0. 18816 116D-04 

0. 18815987D-04 

0. 12937336D-09 

16 

1773344. 

0. 16595599D-04 

0. 16595498D-04 

0. 10072889 D-09 

1? 

1083308. 

0. 14745 156D-04 

0. 14745076D-04 

0. 79625 104D- 10 

18 

1993288. 

0. 131068770—04 

0. 13 186O13D-04 

0 . 63760583D- 10 

19 

2103037. 

0. 11 862329 D-04 

0. 1 1862277 D-04 

0.5 1676748D- 10 
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It should be noted that the value of attraction change of the spher- 
ical Bouguer plate (SBP) at \p * 180° used in Figure 10 was computed by 
using t p = 179.99 degrees. The attraction change of SBP for i p being 
exactly 180° is indeterminate, We can see this by examining equation 
(62) that we get sinhcf> = y which is not defined when i p = 180°. This 
is why we have taken the limiting value when ip * 180°. 

The attraction expressions of the regular and irregular parts of 
the topography are compared by using the rings of the same mass in TDSBP, 
TDPCP and the integration methods. It is found that there is a good 
agreement between the expressions of PCP (equations (65) and (66)) and 
the expression of the irregular part (equation (69) with (68)). Table 
8 is an example for the cylindrical ring 3000 m high. The discrepancy 
between the attraction changes from the two methods is less than 10~ s 
mgal at ip = 1° and it becomes even smaller when ^ increases. 

2.3.3 Conclusion and Discussion 

At the first glance we would say that the expressions for spherical 
Bouguer plate should give a more accurate result than those for the plane 
circular plate since the spherical plate better approximates the earth than 
does the plane plate. However, after examining the numerical results, the 
expressions for the spherical plate may not be appropriate. If we look 
at Figure 9, we see that the potential change at the point P 0 of the spher- 
ical Bouguer plate is increasing while the angular radius of the plate is 
increasing and it would be a positive value for the whole spherical shell 
(^ = 180°). Comparing this to the value from Table 1, the potential change 
at the point P 0 of the spherical shell is negative. Moreover, the at- 
traction change at the point P of a spherical Bouguer plate is not zero 
when ip = 180° (see for example Figure 10). This result contradicts the 
attraction change at the point P of the spherical shell shown in Table 1 . 
This means that the expressions of potentials and attractions that we have 
for the spherical Bouguer plate (equations (36) to (40) and (61) to (64)) 
do not reflect the reality when ip = 180°. So there must be some limitations 
to these expressions. 
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We consider now the expressions for a plane circular plate. As the 
radius of the plane circular plate increases, the value of potential change 
at the point P 0 computed from equations (41) and (43) goes toward the 
potential change of the infinite plane plate. The value of attraction change 
at point P computed from equations (65) and (66) goes toward zero which 
is the value of the attraction change of the infinite plane plate. 

Although the potential change due to condensation at the point P 0 of 
the infinite plane plate is not significantly different from that of the 
plane circular plate, the attraction change at point P on the earth's 
surface of the infinite plane plate is. For example, the difference in 
attraction change between the two models is 0.9 mgals f or the plate of 3000 m 
height and ip - 5° . This magnitude of the attraction change difference 
gives significant effect on the geoid computation. 

When we compared the expressions for potentials and attractions of 
the regular part of topography with those of the irregular part using TDSBP, 
TDPCP and integration methods (explained at the beginning of Section 2.3), 
we found that the expressions of the plane circular plate give a good agree- 
ment both in potential change and attraction change with those of the ir- 
regular part. 

We see that the potential change and attraction change of the plane 
circular plate can be compensated by those of the irregular part. When 
the plane circular plate is extended, its expressions for potentials and 
attractions reflect the correct values. Based on these reasons, the plane 
circular plate will be used to represent the regular part of topography 
in the next chapter. 

In chapter 3, we will develop procedures to compute the potential change 
and the attraction change due to Helmelt's condensation using l°xl°mean 
elevations. 


3 . Calculation of the Indirect Effects and the Attraction Change Effect 
Using l°xl°Mean Elevation Data 

We have seen from equation (11) in Chapter 1 that to compute an accurate 
geoid, we must take into acount (i) the indirect effect, <$N on geoid 
undulation due to the potential change in the second method of Helmert's 
condensation, <5 V discussed in Section 2.1; (it) the effect of attraction 
change in gravity reduction which is equal to the terrain correction, C 
for planar approximation plus a small value of the attraction change of 
the regular part depending on the cap size used in computaion; and probably 
(iii) the secondary indirect effect on gravity, <5 . 

Before we develop a procedure to compute these corrections on the global 
basis, using information of l°xl°mean elevations we should know how geoid 
undulations are practically computed so that, if possible, we can incor- 
porate or apply the corrections properly. 


3.1 Computational Methods for Geoid Undulations 

Stokes' equation in the form of equation (1) cannot be used directly 
in practice because the information on gravity anomalies at every point 
on the earth is not available and the computation on the global basis requires 
a lot of computational effort. Two basic methods are discussed by Rapp 
and Rummel (1975) as follows: 

Method 1 ; the geoid undulation; 

"i ■ N " + Ji? l[ s <*> < A 9l - A W dc, c + "ref (7i > 

where the integral of the second term on the right-hand side is extended 
over cap size a c and the zero order undulation, 

N° = - Ag 0 + Jg-jJ- (75) 
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N ref “ R n " 2 Jo (C ™ C0S " A + 5 ™ sinmX) P nm( sin?) < 76 > 

and - 

Ag ref = y J 2 <n - 1 > J„ < C m cosraX + 5 nm sinmX > P m» (sin?) (77) 

In this method, N re f is the geoid undulation computed from fully 
normalized geopotential coefficients Cfi m and § nm . The short wavelength 
part of the geoid is recovered by using gravity anomalies Ag j with the 
cap size a c . 

Method 2 : the geoid undulation is computed from gravity anomalies in a 

certain cap size a c with Molodensky’s truncation coefficients Q n (tp) 
compensating for the effect of anomalies outside the cap. The basic equa- 
tion is 

Ni = No + // S (ip) Ag j da c + ~ J Q n (^) Ag n (cf>,A) (78) 

a c n= ^ 

where Ag n is the degree anomaly computed from 

n 

Ag n = y(n-l) l (C* cosmX + S nm sinmX) P (sin$) (79) 

m=0 

When equations (74) and (78) are used, the integration over the infin- 
itely small area da c is replaced by the summation over a certain finite 
area, e.g., l°xl°area, then Ag j are mean anomalies over the area. 

In the next section, we consider how to compute indirect effect due 
to potential change from l°xl°mean elevations. 

3.2 Indirect Effect Due to Potential Change 

The total potential change, at point P 0 on the geoid, 6V , caused 
by the second method of Helmert's condensation, can be computed in two parts. 
One belongs to the regular part of topography and the other belongs to the 

irregular part or the topographic variation. Once the cap size \p is chosen, 
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the potential change at the point P 0 of the former part is a straight 
forward computation {see Subsection 2.1.1, page 21). The potential change 
at the point P 0 of the irregular part is computed by summing up V 2 , 

V 3 , etc. terms (See Equations (56) and (49) to (51)). Since these 
V-j terms are in the integral forms which use the continuous function of 
topographic heights, we must modify them to suit the computation which uses 
the discrete mean heights of the block. Now consider the V 2 term from 
equation (49), 

h 3 - hp 

V 2 * - *■ KpR a II da 

0 *0 

Let dS represent the surface area of a sphere of radius R then, 

dS = R 2 do (80) 


quauty 


So we can write 

3 3 

V 2 - - J-kp // dS (81) 

s a; 

In practice, we have to represent the infinitely small area dS by some 
finite area AA-j and we can substitute the integration by a summation. 
Therefore 


V 2 


m 






(82) 


where m = number of AA-j on the surface S 
h.j = average height over an area AA-j 
V-j = distance from the computation point P 0 to the center of the 
small area AA- . 


In the same way, the other terms from equation (50) and (51) can be written 
as 


V 3 = 
V 4 


3 , 5 h5 i - 

a- k p .1. , " 5 - p - " AA i 

1=1 Joi 

. . k p 5 . 

3 3 6 L 


i =1 Aoi 


AA. 
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Up to this point we separate the computation of the potential change 
of the irregular part of topography from that of the regular part. However, 
we have seen from Chapter 2 that if the regular part is represented by a 
plane circular plate, we can then get the same potential change of the plate 
by using the expressions for the potential change of the irregular part 
of topography. This means that instead of using the expression for a plane 
circular plate: (i) we can break the plane circular plate into pieces 

such as consecutive series of cylindrical rings; (ii) use the expression 
for irregular part to compute the effect of an individual piece; (iii) 
sum them up and we still get the same result. 

We are now going to show that by applying the above idea with l°xi° 

mean height information, we can combine the potential change of the regular 
part with that of the irregular part. Instead of breaking the topographic 
masses into pieces of cylindrical rings, we break them into prisms of masses 
above l°xl° areas, i.e., the topography of the earth consists of a series 
of rectangular columns. Now we can use equation (56) with equations (82) 
and (83) to compute the potential change of the plane circular plate where 
AA-j is the area of a l°xl° block. Taking the height of topography h^ 

in equations (82) and (83) to be hp and the height of the reference surface, 

hp in those equations to be zero, we obtain the components of the potential 
change of the -plane circular plate to be: 



We, now, combine the potential change of the regular part (equation (84)) 
and the irregular part (equation (82), (83)) algebraically. Then the indi- 
vidual component of the potential change can be computed from 

m hi 3 

Va= - f kp I (■£“)* AA. 

i = 1 *°i 1 
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(85) 




v ** A* X <^-> 5,AA i 
■ - Hi * l <£>’ • «i 

That is the potential change at the point P 0 on the geoid of the topography 


6V = V 2 + V 3 + V 4 + . , . 

m hf 


= - kp AA i f *- ( p— > ! 
■ ■ ■> 


4 0 \l0j 1 


+ is-(JlL)’ 

'13C4 0 . 1 


(86) 


We can see that the potential change affected by block i depends on two 
parameters h.j and &oj . This makes it easy to compute the potential 
change of the individual block. After summing up the effects of all blocks, 
the indirect effect due to potential change is found from Brun’s formula 
of equation (4). This indirect effect represents the mean value at the 
center of the central block. 


In this section, we have not addressed the problem of the central block, 
which is the block containing the computation point, where io is zero. 

This will be discussed in Section 3.5. Next, we consider how to compute 
the attraction change effect on geoid undulation from l°x 1° mean heights. 


3.3 Attraction Change Effect on Geoid Undulation 

An attraction change effect on geoid undulation is an effect coming 
from the attraction change when the topographic masses are condensed in 
Helmert's condensation reduction. The total attraction change is the sum 
of the attraction changes of the regular and the irregular parts of the 
topography, i.e., 6A = 5A 1 + 6A" . The attraction change of the regular 
part, 6A' = A' - Aj is computed using equation (65) and (66), while the 
attraction change of the irregular part is equal to the terrain correction 
in the planar approximation. Similar to expressions for V 2 , V 3 , V 4 , 
etc. in equations (82) and (83), the terrain correction can be written as 
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9 

i 


<5A 


■■hcp i. M 2 


1=1 


AAj 


Unlike the effect of the potential changes, we cannot combine the effects of 
the attraction change of the regular part and the irregular part into one- 
step computation since the terrain correction is always positive whether or 
not the mass is above or below the surface of the plate. That is the more 
rugged the terrain is, the greater the attraction change is. The sum of 
SA' and iSA" is assumed to be the average attraction change in the central 
l°x 1° block. 

After the total attraction change has been computed, we can add it 
to the free air anomaly before using in Stokes' equation (see equation (9). 
Otherwise, we can compute the correction due to attraction change separately 
(see equation (10)) if the geoid undulation was already computed with uncor- 
rected anomalies. 


In order to compute the attraction change effect on geoid undulation 
separately, we must know how the geoid undulations in the area of interest 
were obtained. For example, if Method 1 of Section 3.1 is used to compute 
the geoid undulation from free-air anomalies and the geopotential coef- 
ficients derived from satellite altimeter. The satellite derived undulations 
are not affected by the topographic masses above the geoid. This means 
that no topographic corrections are needed for the information outside 
cap size c c . Therefore, the cap size used for computing the attraction 
change effect on geoid undulation should be the same as that used in the 
computation of geoid undulation. There are two ways which we can use to 
compute the attraction change effect on geoid undulation from the l°xl° 
mean attraction change that we have. 

Method A ; Use the direct Stokes' equation for the attraction change 
within cap size a c . The effect of the attraction change on geoid undulation 
is then: 

ANi = // 6A S(i|>) da (88) 

a c 

Method B : Use l°x 1° attraction changes over the whole world to devel- 

op correction terms to potential coefficients of the attraction change 
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field using the programs for harmonic development by Colombo (1981), (or 
programs of a similar nature). These correction terms are then used to 
compute the undulation change. The equation for computing the undulation 
change within the cap size a c is modified from equation (51) of Rapp (1981) 
to include the zeroth- and the first-degree components of the attraction 
changes. This equation is: 

ANl <*-0n«"-l»«*n + lr n i 0 l >n 4A n < 89 > 

where 

<$A = l ( 6C hm cosniA + 6S nm sinmX ) p nm< sin ^) (90) 

m=0 

6S " = Jo (SC “ 9) C0Srt + ss ra >SinmX) P rJ si "*> 

6C nrn and 6S nm in equation (90) are corrections in terms of the potential 
coefficients. and are corrections in terms of the anomaly 

coefficients. These coefficients are developed from the attraction changes 
due to the second method of Helmert's condensation. Q n is the Molodensky's 
truncation coefficient of degree n which depends on cap size o c , h 
is the maximum degree of harmonic coefficients used for the computation. 

The derivation of equation (89) is shown in Appendix A. 


3.4 Secondary Indirect Effect 

After we have obtained the indirect effect due to the potential change 
in the l°xl° block, we can use equation (5) to compute the secondary indirect 
effect on gravity, 6 . The value 6 is thus the mean value within the 
l°xl° block too. We can add 6 to Ag + 6A so that we compute the corrected 
geoid undulation. Otherwise, we can compute the correction term due to the 
secondary indirect effect on gravity separately in the same way we compute 
the attraction change effect on geoid undulation (see Section 3.3). 


3.5 Problem of the Central Block 


The central block is defined to be the smallest elevation block which 
contains the computation point P on the earth's surface or point Po on 
the geoid. From equations (86) and (87), we see that <$V" and 6A" vary 
with the inverses of aJj , loj , etc. This means that SV" and <$A" 

are influenced significantly by the information in the innermost zone around 
the computation point Po or P and they vanish when Zoj is very large. 
Consequently we should use more detailed information in the innermost zone. 

For example, if we are working with l°xl° information then we might use 
5'x5' mean elevations in the innermost zone. The size of the central block 
is then 5 1 x 5 ' . We call the l°xl° block which contains the computation 
points as the central zone. The other l 0 xl° blocks within the cap size 
for computing the potential change and the attraction change are called 
blocks in the outer zone. The detailed information however helps to give 
a better value in computing <$A" only. The reason is that in computing 
6V" , if the value of l°xl°mean elevation that we have represents the average 
elevation of all 5 ' x 5 ’ mean elevations in the same l°xl° block, then both 
SV" computed from either size of the mean elevation should be the same. 

It is not necessary then to use 5 ' x 5 ' mean elevations for computing the 
potential change 6V" . The attraction change 5A" of the central zone 
(l°x 1° block) .is zero if only one l°xl°mean elevation is used, The change 

SA" is not zero however if 5 ' x 5 ' mean elevations in the central zone are 

used and at least one of them is different from the l°xl°mean elevation 
of the corresponding block. 

No matter what size of block are used in the calculation of 6V" 

and <$A" in the central zone, we have a problem with the central block 

or zone where z 0 is zero and the effect of the topography is indeterminable. 
Therefore, we now consider how to handle the information of the central 
block. The potential change and the attraction change of the central block 
are computed as follows: 

(a) To compute the potential change 6V" , the central block is sub- 
stituted by the plane circular plate having the same mass as the prism of 
masses above the central block. The following equation is used to obtain 
the radius of the equivalent plane circular plate: 
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Table 9: Postulated l°x 1° Mean Terrain Corrections 


As a Function of the Block Elevation 


Model 1 

Model 2 

Model 3 

Height* 

Ter. Cor. 

Height* 

Ter. Cor. 

Height* 

Ter. Cor. 

0 m 

0 mgals 

200 m 

0 mgals 

300 m 

0 mgals 

200 

1 

400 

1 

500 

1 

400 

3 

600 

2 

800 

2 

600 

5 

800 

3 

1000 

3 

800 

6 

1000 

4 

1500 

4 

1000 

7 

1200 

5 

2000 

5 

1500 

8 

1400 

6 

2500 

6 

2000 

9 

1600 

7 

3000 

7 

2500 

10 

1800 

8 

3500 

8 

3000 

12 

| 

9 

4000 

9 

3500 

14 


10 

9000 

10 

4000 

15 

2400 

11 



4500 

17 

2600 

12 



5000 

18 

2800 

13 



9000 

20 

3000 

14 




< 

3200 

15 





3400 

16 





3600 

17 





3800 

18 





4000 

20 





9000 

21 




*The terrain correction is given for an elevation between the specified 
elevation and the previous elevation. (For example, in model 1, if a 
l°x 1° mean elevation fell between 1500 and 2000 m, the adopted terrain 
correction was 9 mgals). The terrain correction was zero for all 
oceanic blocks. (From Rapp, 1977). 














ORIGINAL PAGE W 
OF POOR QUALITY 


TS (slnfo - slnifi)] 1 

* - 2 s1n [ J 

where 0 is the side of the equiangular block which is bounded by parallels 
i and $ 2 , Both »;/ and 0 are in radians. Note that (sinfo - sin d>i) 
can be approximately computed by cos<hsine, 

(b) The computation of the attraction change or the terrain correction 
of the central block requires a lot of data and efforts to obtain an ac** 
curate result. It can be computed by using detailed height Information 
in a very small block (Dimitri jevich et al., 1976) or by least square col- 
location (Forsberg and Tscherning, 1981). However, the compuation of the 
terrain correction is not the main Interest of this study, so we will use 
the postualted model given by Rapp (1977) to get the mean terrain cor- 
rection of the l°xl° block. The postulated models are shown in Table 9. 

The mean terrain correction given by these models include not only the effect 
of the topography in the central l°xl° block but also the effect of the 
topography in the outer zone. 

3.6 Summary of the Procedures and Equations for Computing an Accurate 

Geoid Using l°xl°Mean Elevations 

In this section we summarize the procedures and the important equations 
used for computing an accurate geoid discussed earlier in the previous sec- 
tions. 

The theoretically more accurate geoid is obtained if Helmert's anomalies 
are used in the Stokes' equation and the indirect effects are taken into 
account. The equation to compute a more accurate geoid undulation is: 

N = Ni + ANj + 6 N 2 + 6 N (11) 

Ni is the free-air geoid undulation. The terms ANi , 6 N 2 and 

6 N are correction terms which can be calculated as follows: 

(a) The potential change effect or the indirect effect, SN : 
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The potential change, <SV is the difference between the potential of the 
topography and the potential of the condensed layer evaluated at point Po 
on the geoid. The effects of the regular and the irregular parts of the 
topography are combined into a one-step computation. A cap size ih is 
chosen for computing the potential change in such a way that the effect 
of ’•he topographic masses outside this cap will be negligible. The potential 
change is the siv of the following two components. 

1, The potential change of the central block is computed using: 


(41) 

(43) 


(42) 

6 is the side of the equiangular block which equals one degree in our case. 


<5V - V' - V) 


where 


hp+Ei 


V - 'rrkp [~hp + h p El + aHn — — ]; E x = 

Vg = 47ikp h R sin |- 


/a z +h^ 

r 


ip = 2 si 


n" 1 r e(sin<p g - slndnl l 

L 4„ J 


and 


a = 2 R sin ^ 


2. The potential change of the blocks in the outer zone is computed from 


m 


5V = -kp l 
i=l 


AA. { 


(in 


- -M— ) 

4 0'Jloj ' 


.} 


(b) The attraction change effect on geoid undulation, AN X . 
We can use either the Stokes' equation: 

ANi = // <$A S(i|») da 

a c 

or the harmonic series: 


(86) 


( 88 ) 


AN, = 


7 l 2 ^ +2 ^I o Q nfi S n 


- 63 - 


(89) 



ORl&WVL PA&E liS 
OF POOR QUALITY 


where n 

W " * Jo (6C ™ cosmX + SS ™ sinmX) P ™ (sin4,) 

6S n ■ j Q < 5C ™’ cosmX + «^ 9 ) sinmX) Pjsln*) (90) 

The attraction change is the difference between the vertical attraction of 
the topography and the vertical attraction of the condensed layer evaluated 
at point P on the earth's surface* A cap size \pz is selected for com- 
puting the attraction change. It may be the same as the cap th depending 
on how accurate the value of the change required. Using a larger cap size 
yields a more accurate result since information in the larger area is taken 
into account. The attraction change consists of the effects of the regular 
and irregular parts of the topography. 


1. The attraction change of the regular part of the topography within the 
cap size \jj z (or the corresponding linear radius a ) is computed from: 

6A ' = A' - Ag 

where 

A 1 = 2tt k p (a + h p - /a 2 +hp ) (65) 


A! = 27Tkph„ (1 2 — ) < 66 ) 

s p /a^TFP - 

2. The attraction change of the irregular part of the topography or the 
terrain correction is computed from: 


<$A" = 


ikp 


m 

i 

1=1 


( h j - hp ) 2 


£ 3 o 


AA. 

l 


(87) 


Rapp's postulated mo * 4 el 1 given in Table 9, however, will be used to obtain 
the terrain correction in this study. 


( c ) The secondary indirect effect on qeoid undulation, 6 N 2 • 


where the secondary indirect effect on gravity, 6 is obtained from 

6 = 0.3086 <SN (5) 

In the next chapter, we present the results of indirect ; ‘fects due 
to the potential changes, attraction change effects on geoid undulation 
and the secondary indirect effects in the test area using procedures ex- 
plained in this section. 


4. Results in the Test Area 

An area referred to as an output area, which is used to test the previous 
ideas, was selected to extend from latitude 20°N to latitude 50°N and from 
longitude I30°W to 65°W. A data area where mean elevations of l°x 1° equi- 
angular blocks were used to compute the potential and attraction changes, 
then, has to be larger than 30° x 65°. It depends separately on the cap sizes 
sizes that we use to compute the changes and the cap size o c in Stokes' 
equation. For example, Figure 11 is the diagram showing the data area cf 
l°x 1° mean elevations needed to get attraction change effects on geoid 
undulation in the output area of 30° x 65° when the 10° cap is used in Stokes' 
equation and the 5° cap is used to compute the attraction changes. The effect 
of the potential change on the geoid undulation is computed directly from the 
the indirect effects due to potential change within the same output area 
is smaller. Figure 12 is the contour map showing topographic heights above 
the geoid in the. output area. 

In the next section, we look at the effects of the potential changes 
due to Helmert's condensation in the test area. 


4.1 Effect of Potential Changes on Geoid Undulations 

The indirect effect due to the potential change in the second method 
of Helmert's Condensation was computed by the procedure explained in section 
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geoid have to be condensed onto the geoid in the second method of Helmert' 
condensation. The contour interval of the map is 200 meters. 




3.6(a), It was anticipated that the topographic masses further away from the 
computation point would have a very small contribution on the total potential 
change and thus the indirect effect. To test how the remote topographic 
masses have an influence on the indirect effect, cap sizes of 3°, 5° and 
10° were used to compute the indirect effects at the center of each l°xl° 
block in the test area and compared with those when the 20° cap was used. 

It was found that the differences are in the order of few tenths of a millimeter. 
The indirect effect computed using the masses only in the central l°xl° 
block (called the 0° cap) were also used in the comparison. The maximum 
discrepancy was only 4 millimeters. The statistics of the indirect effects 
in the test area for the cap sizes of 0°, 5° and 20° are shown in Table 
10. We can see that the remote topographic masses have insignificant con- 
tributions on the indirect effect due to the potential change at the centi- 
meter level of accuracy. The insignificant contributions of the remote 
masses are understandable by examining the V 2 term of equation 
( 86 ). The V 2 term consists of the factor (- I 11 -) 3 . If the mean elevation 

~ 0 i' j"| j 

of the block next to the central block is 1000 meters, then — of this 
block is less than 10" 2 . By neglecting the mass of this block, the relative 
error in computing the indirect effect is less than 1CT 6 . Presumably, 
masses in the block further away from the central block cause a smaller 
relative error. As a result, the accumulation of relative errors in neglecting 
all masses around the central block is still very small and the effect of 
these masses is insignificant. From this investigation, we can say that 
only the mean elevation of the central block (block containing the compu- 
tation point) is needed for computing the indirect effect due to the poten- 
tial change. 

To illustrate the indirect effect in the test area, the mean value 
at the center of each l°xl° block was used to produce a contour map. The 
contour map of the indirect effects is shown in Figure 13. The maximum 
absolute value of the indirect effect in the test area is 51 cm. It occurs 
at the highest point in the area. We can observe the high correlation between 
the topography and the indirect effects by comparing Figure 12 and Figure 
13. We see that the peaks of the indirect effects occur at the same locations 
where the peaks of topography occur. We can see also that the indirect 
effects due to potential changes in the ocean area are zero. 
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Table 10 : The indirect effects (6 V/y) in the test area f 

using different cap sizes (20°, 5° and 0° for 
the calculation of the potential changes. 


Statistics 

— 

Capsize 

20° 

5° 

0° 

mean (m) 

-0,0343 

-0.0342 

-0.0340 

S.D (m) 

0.0728 

0.0728 

0.0721 

max (m) 

0.0000 

0.0000 

0.0000 

min (m) 

-0.5076 

-0.5076 

-0.5036 

RMS (in) 

0.0805 

0.0805 

0.0798 


4.2 Attraction Changes and Their Effects on Geoid Undulation 

The attraction change due to the second method of Helmert's condensation 
for 1^x1° block was computed by the procedure explained in Section 3.6(b). 

In order to study the effects of remote masses on the attraction changes, 
cap sizes of 10° and 5° were used in the computation and compared with 
those when the 20° cap was used. Rapp's postulated model No. 1 for esti- 
mating the attraction changeof the irregular part of the topography or the 
terrain correction was used in this study. Comparison of the three cap 
sizes is shown by the statistics in Table 11. The values under column "Dif- 
ference between cap sizes" were obtained by subtracting the attraction change 
computed using one cap size from the other of the same l°x 1° block, then 
the statistics of the differences in the area were computed. From the table, 
we see the maximum absolute value of the difference from the 20° cap is 
0.24 for the 10° cap and 0.69 for the 5° cap. Based on the most accurate 
values for the attraction changes that we obtain from the 20° cap, the root 
mean squares (RMS) of the differences of the 10° cap and 5° cap indicate no 
significant difference. In the mountainous terrain we have the larger 
magnitude of difference between the different cap sizes than that in the flat 
area. We must decide what order of magnitude differences can be tolerated 
in our work. In this study, the 10° cap was used for further investigation. 
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The contour interval of the map is 2 centimeters. 
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Table 11 ; Comparison of the attraction changes in the test area 
using different cap sizes (20°, 10° and 5°) for the 
computation. 


Statistics 

Cap Sizes 

Difference between 
cap sizes 

20° 

10° 

5° 

20°and 10° 

20° and 5° 

mean (mgal ) 

2.1006 

2.0848 

2.0551 

-0.0158 

-0.0455 

S.D, (mgal) 

3.0510 

3.0223 

2.9685 

0.0332 

0.0964 

max. (mgal) 

11.7394 

11.5009 

11,0517 

0.0000 

0.0000 

min, (mgal) 

0.0013 

-0.0025 

-0.0047 

-0.2386 

-0.6877 

RMS (mgal) 

3.7036 

3.6710 

3.6098 

0.0370 

0.1067 


Figure 14 is the contour map of the attraction changes in test area 
due to the Helmert's condensation method. This contour map was produced 
from l°xl° grid values at the center of the blocks. Since the total attrac- 
tion change is dominated by the terrain correction which is directly dependent 
on the elevation of the computation point, it can be anticipated that the 
peaks of the attraction changes occur where the peaks of the topography 
do. From this map, the areas having steep slopes are accentuated by the 
high gradients' of the attraction changes. Similar to the potential change, 
there are no attraction change in the ocean area. 

The attraction changes obtained from the previous step were used to 
compute the attraction change effect on geoid undulation. First we used 
Stokes' equation with cap size, a c of 10° (Method A in Section 3.3). There 
were two parts of the attraction change being computed, i.e., the effect 
of the regular part and the irregular part of the topography. Since the 
regular part has a small contribution on the total attraction change, we 
investigated how significant the attraction change of the regular part is 
in computing the correction to the geoid undulation. Two models of the 
attraction change were then used for this purpose. Model 1 consists of 
the attraction changes of both the regular and the irregular parts of the 
topography. Model 2 contains only the attraction change of the irregular 
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Figure 14 : The contour map of the attraction changes due to the second method of Helmert's 

condensation. The attraction chanoe was computed at the center of each l°xl° 
block using Rapp's postualted model no. 1 (Table 9) and equations ( 65) and ( 66). 
These attraction changes were then used to produce the contour map. The contour 
interval of the map is 0.5 mgal. 
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part or the terrain correction. The result of using Stokes' equation (88) 
to compute the attraction change effect for geoid undulation is shown in 
Table 12. The maximum value in the area is 9.50 m computed by using model 
1. Dimitrijevich et al. (1976) use an average value of the terrain cor- 
rection over some defined cap to compute its effect on height anomalies. 
Their result for the 10° cap is the following: 

Aver, terrain 

correction 1 2 3 4 5 6 7 8 9 10 12 15 

(mgal ) 

Height 

anomaly 1.3 2.7 4.0 5.4 6.7 8.1 9.4 10.8 12.1 13.4 16.1 20.2 

(m) 

Considering this and the attraction changes shown in Figure 14, the maximum 
value of 9.50 m we have is not an unexpected magnitude which can be obtained 
for the 10° cap. 


Table 12 : Comparison of the attraction change effects on geoid undulation 
computed by using different models of the attraction change 
in the Stokes' equation (Method A) with a = 10°. 


Statistics 

Corrections to geoid undulation 
due to attraction change 
(Method A equation (88)) 

Difference 

between 

Models 

Al & A2 

Model 1 

Model 2 

Mean (m) 

2.7523 

2.8361 

0.0838 

S.D. (m) 

2.7887 

2.8938 

0.1077 

Max. (m) 

9.4970 

9.9111 

0.4318 

Min. (m) 

0.0000 

0.0000 

0.0000 

RMS (m) 

3.9176 

4.0513 

. .. . 

0.1365 


Al in Table 12 refers to the procedure whose attraction change effects 
on geoid undulation are computed by using model 1 of the attraction change 
in method A. If we take the difference between the undulation changes computed 
by method A using model 1 and model 2, and compute the statistics, we obtain 
the values in column "Difference between models Al and A2." It can be seen 
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using the attraction changes of both the regular and irregular parts of the 
topography (model 1) in the Stokes' equation (88) with a = 10° - The l°x 1° 
grid values of the attraction change effects were then used to produce the contour 
map. The contour interval of the map is 0.5 meters. 
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that RMS of the differences is 13.6 cm when we use model 2 instead of the 
more accurate model 1. Therefore, model 2 cannot be used as an approximation 
of attraction changes to compute the attraction change effects on geold 
undulation for the geold of 10 cm accuracy level. The contour maps of the 
corrections to the geoid undulation due to the attraction changes computed 
by Stokes’ equation (88) using models 1 and 2 are shown in Figures 15 and 
16 respectively. Comparing these two figures, we can see the major difference 
only in the highest contour of 9.5 m. The other contours are almost the 
same. We can see that the effect of the big mass in the high mountain area 
on the corrections to the geoid undulation can extend over the ocean area 
even if there are no attraction changes in the ocean area. 

The l°x 1° attraction over the whole world, were also used to develop 
a set of spherical harmonic coefficients. The set of spherical harmonic 
coefficients was developed into potential coefficients and anomaly coefficients 
(if n < 2) up to degree 180 (Colombo, 1981), These coefficients were then 
used to compute the corrections to the geoid undulations in the test area 
using equation (89). The correction computed by method B is shown in Table 
13. 


Table 13 : Comparison of the corrections in gravity reduction computed 
from the harmonic coefficients (method B) with cr ~ 10° and 
comparison between method A (see Table 12) and method B. 


Statistics 

AN by method B 

Differences of AN 
between methods 

Model 1 

Model 2 

Difference 

A1 & B1 

A2 & B2 

Mean (m) 

2.5469 

2.9987 

0.4519 

-0.2055 

0.1626 

S.D. (m) 

3.1513 

2.9440 

0.3881 

0.5230 

0.2944 

Max. (m) 

9.8689 

10.1256 

1.8309 

0.6204 

0.7758 

Min. (m) 

-1.9339 

-0.1768 

0.0934 

-1.9339 

-0.5501 

RMS (m) 

4.0512 

4.2018 

0.6001 

. .... 

0.5617 

0.3363 
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The contour map of the attraction change effects on geoid undulation, pro- 
duced from l°x 1° grid values which were computed by method B via the spher- 
ical harmonic series (equation (89)) using model 1, isshownin Figure 17. 
Comparing Figure 17 to Figure 15, the patterns of the contours are the same. 
There are some contours extending over the ocean area on the east side of 
the continent, 

Also shown in Table 13 are the differences of the corrections to geoid 
undulation due to the attraction change computed by methods A and B when 
the same model of the attraction change was used. Figure 18 is the contour 
map showing the differences of the corrections when model 1 of the attrac- 
tion change was used. The differences between the two methods are system- 
atic. The differences are caused by omitting the harmonic coefficients 
of degrees higher than 180. It was found that after the l°x 1° attraction 
changes were transformed into the set of spherical harmonic coefficients 
(potential coefficients and anomalies coefficients (if n < 2)) up to degree 
180, the attraction field was already distorted. The distortion of the 
field was detected by recreating the attraction changes from the set of 
the harmonic coefficients and comparing with the original attraction changes. 
Table 14 is the comparison by statistics in the test area between the orig- 
inal and the recreated attraction changes of model 1. The RMS difference 
is 0.6 mgal and the maximum difference is almost 3 mgal which is quite a 
large value. 

Table 14 : Comparison of the original attraction changes (model 1) and 

the recreated ones computed from the sec of the spherical harmonic 
coefficients up to degree 180, (potential and anomaly coefficients 
(if n < 2)). 


Statistics 

Original Attrac- 
tion changes 

Recreated Attrac- 
tion changes 

Differences 

Mean (mgal) 

2.0975 

1.9481 

-0.1493 

S.D. (mgal) 

3.0353 

3.1801 

0.5951 

Max. (mgal) 

11.5884 

10.4839 

2.0466 

Min. (mgal) 

-0.0017 

-1.6459 

-2.8877 

RMS (mgal) 

3.6888 

3.7287 

0.6134 
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We furthered investigated how well equation (89) agrees with the Stokes' 
equation (88). The recreated attraction changes were therefore used to 
compute the corrections to the undulations due to the attraction changes 
by equation (88). The results were then compared with the corrections computed 
by equation (89). The statistics of the comparison in the test area is 
shown in Table 15. The mean difference is only 4.5 cm. This shows a good 
agreement of equation (89) and the Stokes' equation (88). > 


Table 15 : Comparison of the attraction change effects on geoid undulation 
using equation (89) and the Stokes' equation with the recreated 
attraction changes from the harmonic coefficients. 


Statistics 

AN by 

Stokes' Eq. (88) 
with the recreated 
attraction changes 

AN by 

harmonic series 
(equation (89)) 

Differences 

Mean (m) 

2.5922 

2.5469 

0.0453 

S.D. (m) 

3.1330 

3.1513 

0.1418 

Max. (m) 

9.8583 

9.8689 

0.4943 

Min. (m) 

-1.8844 

-1.9339 

-0.4797 

RMS (m) 

4.0657 

4.0512 

0.1488 


4.3 Secondary Indirect Effect 

From section 4.1 we have seen that the maximum and the mean indirect 
effects on geoid undulation, 6N in the test area are 0.51 and 0.03 meters 
respectively. If we use equation (5) to compute the secondary indirect 
effect on gravity 6 , we will have the maximum effect of only 

0.15 mgal and the mean of 0.01 mgal . The magnitude of attractions in this 
order does not give any significant effect on the value of geoid undulation 
at all. Therefore, the computation of the secondary indirect effects on 
geoid undulation was neglected. 
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topography (model 1) in the spherical harmonic series (equation 89) with the 
10° cap. The l n xl° grid values of the attraction change effects were then 
used tG produce the contour map. The contour interval of the map is 0.5 meters 








Figure 18 : The contour map showing the differences of the attraction change effects 
on geoid undulation computed by method A (equation (88)) and method B 
(equation (89)). Model 1 of the attraction change was used in computing 
the attraction change effect. The contour map was produced from l 0 xl°grid 
values. The contour interval of the map is 0.1 meters. 




4.4 Comparisons between Doppler and Topographically Corrected Gravimetric 
Undulations 


I 
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The calculations derived in the previous sections represent theoretical 
concepts that should be subject to some actual verification. One way to 
do this is to compare geoid undulations implied by Doppler satellite obser- 
vations at selected sites with those gravimetrically computed undulations 
based on data with and without the correction described in this study. A 
computation of the undulations, without correction terms, has been described 
by Lachapelle (1979) at 65 stations in the US using the GEM10B potential 
coefficients and mean surface gravity anomalies {uncorrected for terrain 
effect) of varying sizes of anomaly blocks i.e., 5 ' x 5 ' , 30 1 x 30 1 and l p xl° 
blocks. These undulations are referred to as uncorrected gravimetric 
undulations. 

Using the procedures of this paper, the correction terms were evaluated 
at each of the 65 stations. The correction terms applied to the station 
are the mean values of the block where the station is located. The attrac- 
tion change corrections were computed using the direct integration procedure 
(equation (88)) and the spherical harmonic technique (equation (89)) based 
on l°xl°mean attraction changes (the same values obtained in Section 
4.2). By using. l°x 1° attraction changes for computing the corrections, 
an inconsistency in the computation is created due to the fact that Lachapelle 
(ibid.) used 3 sizes of anomaly blocks. The inconsistency can be remedied 
if we were to use 5'x5' mean elevations to compute the attraction changes. 


The corrections were used to compute a theoretically more correct 
undulation which is referred to as a topographic corrected gravimetric 
undulation. The corrected undulation was then compared with the Doppler 
derived undulation given by Lachapelle (ibid.). Before the actual 
comparisons are made, it is critical to make sure the Doppler undulations 
and the gravimetric undulations refer to the same ellipsoid. 


In order to do this we first express the Lachapelle gravimetric undu- 
lations in the form: 

N g,L * fi o + R g 
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where R. is the zero-order undulation of the geoid so that N_ . will f 

0 9 » L 

refer to a specific ellipsoid. Values of R g were computed by Lachapelle 

using constants of the Geodetic Reference System 1967 (GRS67) . Lachapelle 

computed R g and applied a correction of +25m (6378160-6378135m) 

to make the gravimetric undulations compatible with the ellipsoid used 

in the determination of Doppler undulations. Lachapelle then computed 

the mean difference aN defined as: 


n 



where N Q is the Doppler derived undulation. Values of AN were determined 
using several sets of potential coefficients and anomaly sets. For our 
purpose we use the computations with the combination of GEM10B and the 
anomaly set. In this case the value of /IN was -25.7 meters. Now Lachapelle 
gives the values of R' which are defined as: 

* L + .25m + AN 
s y 

Thus we can re-construct the values of Rg from the above equation: 

R g = Rg - 25m - AN (92) 

Note that the mean values of Rg would be zero with respect to the Doppler 
undulations. 

In order to utilize these R g values we must compute the zero-order 
undulation of the geoid so that undulation given by equation (91) will 
in fact refer to the Geodetic Reference System 1967. From Rummel and 
Rapp (1976, equation (25)) we have: 

N 0 - «a[l + *(i|> c )] - f(t|. c ) 
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where <$a is the difference between the equatorial radius of the mean 
earth ellipsoid and that of GRS67, k<SM is the difference between the 
geocentric gravitational constant of the mean earth ellipsoid and that 
of GRS67, and <i>(if) ) is defined as: 

88 47 H s (^) d a c 

a c 

Based on the recent estimates of the parameters of the mean earth 
ellipsoid, we find (see Appendix C) 

R Q = -26.11 meters 

Thus the geoid undulation with respect to GRS67 is the sum of and 
fl given by equation (92). Thus: 

Ng.GRSe? = - 26 ’ 11 "’ • 25m + N g - 4N 

For comparisons with the Doppler undulations a correction of the above 
undulation is needed so that the values refer to the same ellipsoid. 

We would have 

R g,D = S g,GRS67 + (6378160 ' 6378135m) 

Thus 

N n = -26.11m - 25m + 25m + N ' - AN 
g,D „ g 

= -26.11m + Rg - AN 

For the AN of -25.7m we have 

R n n = R' - 0.41 meters 

y > u g 

We now are able to compare the values of R g D with the undulations 
derived from the Doppler stations. Such comparisons are given in the 


s 

/ 


first numerical column of Table 16. We see a systematic difference (Doppler 
minus gravimetric N ) of 0.4 meters with a difference standard deviation 
of ±2.2 meters. At this point we apply the correction terms based on 
the Helmert's condensation procedure to the uncorrected (for topography) 
undulation given by equation (93). The details of the corrections at 
the65US stations are shown in Appendix B. The results of the comparisons 
are given in the last two columns of Table 16. 


Table 16: The differences: Doppler undulation minus gravimetric 
undulations. The Doppler undulations are those from 
Lachapelle (1979). 


Statistics 
of the 
Differences 

Uncorrected 
Undulations 
Derived from 

Lachapelle (1979) 

. . 

Corrected Undulations 

When AN 
Computed by 
Equation (88) 

When AN 
Computed by 
Equation (89) 


0.41 

-3.42 

-3.36 


2.22 

3.12 

3.25 


5.83 

5.03 

5.23 


-3.78 

-8.91 

-9.29 


2.24 

4.61 

4.66 


From Table 16, we can see the deterioration in the agreement between the 
gravimetric undulations and the Doppler undulations in all aspects, after 
the corrections were added. The absolute mean differences, the difference 
standard deviations and the difference RMS of the corrected sets are larger 
than the corresponding values of the uncorrected set. 

In performing these comparisons we have assumed that the Doppler 
derived coordinates are geocentric since the gravimetrically determined 
values should theoretically be so. In fact there is a good evidence that 
there is a 4m z-axis shift of the Doppler origin with respect to the geocenter. 
Thus we should correct the original Doppler undulations to convert them 
to a true geocentric systems. We have: 

Np = N D + Az s i n d? 
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where az is the z-axis shift of the Doppler origin with respect to the 
geocenter, and <j> is the geodetic latitude of the Doppler station. We 
now repeat our gravimetric-Doppler comparisons using the N^ values. 

The results are given in Table 17. 

Comparing Table 17 to Table 16, we see the decrease in the values 
of the difference standard deviations in all three sets of the gravimetric 
undulations. This indicates that the systematic effect due to z-shift 
of the Doppler origin with respect to the geocenter does exist. 

Table 17: The differences: Doppler undulations minus gravimetric 

undulations. The Doppler undulations are those from 
Lachapelle (1979) corrected for systematic effect due 
to 4m z-axis shift. 


Statistics 
of the 
Differences 

Uncorrected 
Undulations 
Derived from 
Lachapelle (1979) 

Corrected Ur 

idulations 

When AN 
Computed by 
Equation (88) 

When AN 
Computed by 
Equation (89) 

Mean (m) 

2.86 

-0.97 

-0.91 

S.D. (m) 

1.99 

2.87 

3.02 

Max. ( m ) 

7.87 

6.78 

6.98 

Min. (m) 

-0.90 

-6.09 

-6.63 

RMS (m) 

3.48 

3.01 

3.13 


When the absolute mean differences of the corrected undulations are 
compared to that of the uncorrected undulations in Table 17, the improvement 
due to the topographic corrections is obvious. The standard deviations 
of the two corrected sets are however larger than that of the uncorrected 
set. The poorer standard deviations may be caused by the inaccurate values 
of the terrain corrections used in the computation of the attraction change 
effect on geoid undulation. 

5. Summary and Conclusion 

The indirect effects on the geoid computation due to the second method 
of Helmert'r condensation were studied. Theoretically, the use of Helmert's 
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anomalies in Stokes' equation should give a more accurate geoid than the 
use of uncorrected free-air anomalies. Using Helmert's anomalies gives 
three types of corrections to the free-air geoid. First, the indirect 
effects on the geoid undulation are the result of the potential changes 
in condensing the topographic masses above the geoid onto the geoid in 
Helmert's reduction method. The second one is the correction to geoid 
undulation affected by the attraction changes in the reduction process. The 
third one is the secondary indirect effects on geoid undulation which is the 
consequence of the indirect effects caused by the potential changes. 

The computation or evaluation point for computing a potential change 
in a gravity reduction is not at the same location as the computation point 
for computing an attraction change. In the second method of Helmert's con- 
densation, the potential change is the difference between the potential 
of the topography and the potential of the condensed layer, evaluated at 
some point P Q on the geoid . The attraction change, on the other hand, 
is the difference between the vertical attraction of the topography and 
the vertical attraction of the condensed layer. It is evaluated at a point 
P on the earth's surface . 

In the derivation of potential changes and attraction changes, the 
topographic masses above the geoid are divided into a "regular" part and 
an "irregular" part (Moritz, 1968). The "regular" part of the topography 
is sometimes regarded as an approximation to the actual topography. Three 
models can be used to represent the regular part of topography: the spherical 

Bouguer plate, the plane circular plate and the infinite plane plate. The 
"irregular" part of the topography consists of the topographic masses which 
deviate from the plate surface. 

The first part of this investigation was an investigation of the poten- 
tial change and the attraction change of the regular part of topography as 
computed from the three different models. These changes were then compared 
with the irregular part when the same amount of mass was used. We found 
the closest agreement between the results from the plane circular plate 
and the irregular part of the topography. 
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Next we developed computation techniques to compute the potential 
changes and the attraction changes from l°x 1° mean elevations. A test area 
covering the United States was then used to investigate the actual effects 
of the three corrections to the free-air geoid. It was found that the maximum 
indirect effect on geoid undulation (potential change effect) in the test 
area was 51 centimeters. The maximum attraction change effect on geoid un- 
dulation was 9.50 meters when the 10° cap was used in Stokes' equation (88). 

The attraction change effect on geoid undulation computed from the set of 
harmonic coefficients (equation (89)) up to degree 180 was not the same as 
that computed by Stokes' equation. The discrepancy between the two methods 
was mainly caused by the higher degree coefficients being omitted, We have 
shown that equation (89) produced results approximately equal to Stokes' 
equation (88) when the same attraction change field was used in both equations. 

The secondary indirect effects on gravity were of small magnitudes in 
the test area* the mean value in the area being 0.01 mgal . Thus, the effects 
could be neglected in the computation. 

Doppler and gravimetric undulations reported by Lachapelle (1979) were 
used to investigate how well the corrections due to the potential changes 
and attraction changes could improve the gravimetric undulations. The 
gravimetric und,ulat-ions were slightly modified so that they would be comparable 
to the Doppler derived undulations, i.e., they would refer to the same 
ellipsoid (section 4.4). These undulations are referred to as the uncorrected 
gravimetric undulations. The corrections were then applied to the uncorrected 
undulations to obtain the so called topographically corrected gravimetric 
undulations. 

Since there is a tendency that the Doppler origin may have a z-axis 
shift of 4 meters with respect to the geocenter. The corrections due to 
this shift were computed and applied to Lachapelle's Doppler undulations 
to get the "corrected" Doppler undulations. The uncorrected and the 
corrected gravimetric undulations were then compared with the original and 
the corrected Doppler undulations. 
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Our comparisons confirm that there exists the z-shlft In the Doppler 
origin. Based on this fact, the corrections due to the potential changes 
and the attraction changes in the Helmert's condensation, show that 
systematic differences are reduced although the dispersion of the results 
increases. More specific conclusions depend on the determination of more 
accurate terrain correction values. 

The terrain correction used in this study directly depended on the 
elevation of the evaluation point. The type of terrain surrounding the 
evaluation point was not considered, This assumption leads to a remarkable 
error in the terrain correction when the evaluation point is located on 
a high plateau. In such an area, the terrain correction would be very small, 
or zero, if equation (87) were used. Accurate values for the terrain cor- 
rection may be obtained if mean elevations in smaller blocks were used as 
done in Dimitri jevich et al . (1976) ( 5 ' x 5 ' blocks) and Curiale et al. (1981) 

( 5/8 1 x 5/8' blocks). 

If the more accurate corrections due to the attraction change still cannot 
improve the agreement of the Doppler undulations and the corrected gravi- 
metric undulations, it is suspected that Helmert's anomalies may not be 
appropriate for geoid computation in the high mountain area such as the Rocky 
Mountains. 

It is recommended, then, that the use of Molodensky's correction in 
mountainous regions be studied to see if a better agreement between gravi- 
metric and Doppler undulations may result. Two comments should be made 
in regards to solving the Molodensky boundary value problem. One is that 
instead of the geoid undulation, the height anomaly is obtained from the 
solution. It is stated by Heiskanen and Moritz (1967, p. 294) that these 
two are very close in magnitude. The second point is that the linear de- 
pendency of the gravity anomaly on the elevation must not be presumed. The 
reason for this is that the assumption will lead us to the first-order 
correction (Gj term) which is the same as the terrain correction (Moritz, 

1980, p. 414). 
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Appendix A 

Derivation of Fquation 


ORIGINAL 

01? poOR QUALITY 


Ml 


Two methods for computing the geoid undulation are presented in section 
3.1. The undulations computed by both methods must be the same. If we 
equate equations (74) and (77), we get 

// (Agr - Ag re f) S(^) da c + R A n 

03 

= // Agi sm da c + J Q n Ag n (Al) 

c '~ c 

where 

A n = Jo (C ™ cosmA + 5 ™ s1nm! ' ) P ™> (sin * ) 

Ag n = v(n-i) A n 


The maximum degree n of the harmonic series is usually selected in such 
a way that the error in neglecting the higher degree terms (n > n) will 
be very small. The upper limit of the summation in the second term on the 
right hand side of equation (Al) then can be set to n . The lower limit 
of this summation must be extended to zero because the zeroth- and the first- 
degree components of Ag do exist in the set of the harmonic coefficients 
developed from the attraction changes. Since the zeroth and the first- 
degree potential coefficients are not defined, we have to use the anomaly 
coefficients to compute the components Ag 0 and Agi . Equation (Al) can 
be written as: 


4y // < A 9* " A 9 re f) da c + R [ 0 A 


n=2 


Sv ff A 9* S ^) dCT r + ST l Qn A 9 


Tiy 

w c 

VmJ!^ref SW,) da c 


c ^ n =0 n n 


°c 


« ! Vr; ! I Mn-n A J + 

n=2 n n=2 n 


R l 

J T '6 1 y Q Ag 
n ^y n =o n yn 


■ T X < 2 ' A n + Ty j 0 Q n 
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-3 


(A2) 



ORKartAS. ESS 
OR POOR QUALITY 

If ^ m 9) and 5^) are anomaly coefficients, then we can compute 
9 n by: 

i9n = jo <C ™ 9) C ° 5mX+ 5 ™ 9) sinm *> W sl ' nl) 


Thus equation (A2) can be used to compute the undulation from harmonic 
coefficients within the cap size of a c when the zeroth and the first- 
degree components of Ag exist. 


To transform equation (A2) to equation (89), the following list is 

used: 



Ag 


SA 


A n 


*n 


«n 


sr n 


C m • 
C<Ag) 
Sm ’ 

^nm 

$(Ag) 

nm 

^nm ’ 

«C (A9 > 
nm * 

6S 

nm 

«s (49) 

nm 







Appendix B 

Comparison between Doppler and Topographic Corrected 
Gravimetric Undulations 


V 






9 

i 


Notations used in the table: 

N d = doppler derived undulations (from Lachapelle, 1979) 

r = uncorrected gravimetric undulations (see Section 4.4) 

DN = N D - N g 

6N = correction due to the potential change 

AN 1 = correction due to the attraction changes computed by Stokes' 
equation (88) 

fig 1 = the toTected gravimetric undulation when aN 1 is used 
= fTg + 6N + AN * 

DN 1 - N D - fig 1 

AN 2 = correction due to the attraction changes computed by the harmonic 
series (equation (89)) 

flg 2 =f!g + <5N+AN 2 

DN 2 = N d - iVg 2 
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Appendix C 

The Zero Order Undulation 
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The zero order undulation for the geoid computation by Method 1 and 2 
of section 3.1 is given by (Rummel and Rapp, 1976, equation (25)): 

R 0 = <Sa 1 + 4>(^ c ) - 4>(ifi c ) ( C1 ) 

With 

= .^f J7 s W dcJ c (C2 ^ 

a c 

where 5 a is the difference between the equatorial radius of the mean earth 
ellipsoid and that of the adopted reference ellipsoid, and k 6 M is the differ 
ence between the geocentric gravitational constant of the mean earth ellipsoid 
and that of the adopted ellipsoid. 

Note that in equation (Cl) is not that as N Q given in equation 
(75). This is because the global mean anomaly, Ag Q must be subtracted from 
Ag j in the second term of equations (74) and (78), to assure that the gravity 
anomalies used in Stokes' equation have a global average equal to zero. See 
the detailed derivation in Rummel and Rapp (ibid). 

At the present, the best-known constants for the mean earth ellipsoid 
which will be used in our calculations are: 

(from Moritz, 1979), 

kM = 3.9860047 x 10 14 m 3 /sec 2 
(from Rapp, 1982), 
a = 6378135.6 m. 

The constants i or the Geodetic Reference System 1967 (GRS67) are: 

kM = 3.98603 x 10 14 m 3 /sec 2 
a = 6378160 m. 

That is !<<$M = -253 x IQ 7 m 3 /sec 2 and 5 a = -24.4 meters. 
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For a cap of 10° in Scokes' integral f 4>(^ c = 10°) s 0.2068. Using 
R ® 6371000 meters and y = 9.8 m/sec 2 , the zero order undulation of 6RS67 
is: 

= -26.11 meters. 
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